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ABSTRACT 


A  research  and  development  study  of  stochastic  finite  element  methods  (SFEM)  for  structural 
reliability  analysis  is  presented.  A  detailed  state-of-the-art  literature  review  of  various  SFEM 
methodologies  and  the  diverse  technical  subjects  relating  to  probabilistic  finite  element  analysis  is 
summarized.  Case  studies  of  two  existing  SFEM-based  reliability  analysis  computer  programs  were 
conducted,  with  a  discussion  of  the  methodology,  main  features,  structure,  and  capabilities  of  both 
systems.  The  design  and  implementation  of  a  software  package  for  random  field  discretization 
(RANFLD)  and  SFEM-based  reliability  analysis  (STOVAST)  is  then  presented.  The  package, 
developed  on  the  basis  of  the  finite  element  reliability  method  (FERM)  approach,  and  designed  for 
operation  with  the  VAST  commercial  finite  element  analysis  system,  was  verified  using  several 
problems.  Also,  in  this  work,  the  quantification  and  analysis  of  bias  and  modelling  errors/uncertainties 
in  engineering  analysis  are  discussed.  Sources  of  bias  and  modelling  errors/uncertainties  in  SFEM- 
based  reliability  analysis  are  identified  and  suggestions  are  offered  for  the  application  of  realistic 
formulations  for  accounting  for  this  class  of  uncertainties.  Finally,  recommendations  are  given  on 
requirements  for  further  development  work  that  will  broaden  the  analysis  scope  of  the  STOVAST 
system  and  permit  straightforward  applications  to  ship  structures. 


RESUME 


On  presente  une  etude  de  recherche  et  developpement  portant  sur 
1 ' application  a  la  fiabilite  structurale  des  methodes 
stochastiques  a  elements  finis  (MSEF) .  On  resume  une  analyse 
bibliographique  detaillee  des  connaissances  de  pointe  sur 
diverses  MSEF  et  les  divers  aspects  techniques  relies  a  1 ' analyse 
probabiliste  par  elements  finis.  On  a  effectue  des  etudes  de  cas 
de  deux  programmes  inf ormatiques  d' analyse  de  f iabilite  basee  sur 
les  MSEF  et  on  a  examine  la  methode,  les  caracteristiques 
principales,  la  structure  et  les  capacites  des  deux  systemes .  On 
traite  ensuite  de  la  conception  et  de  la  mise  en  application  d'un 
progiciel  de  discretisation  a  champ  aleatoire  (RANFLD)  et 
d' analyse  de  fiabilite  basee  sur  les  MSEF  (STOVAST) .  Le 
progiciel ^base  sur  la  methode  d' analyse  de  fiabilite  par  elements 
finis  (MFEF)  et  destine  a  etre  utilise  avec  le  systeme  commercial 
d1 analyse  par  elements  finis  VAST  a  ete  verifie  par  la  resolution 
de  plusieurs  problemes.  Ce  travail  traite  egalement  de  la 
quantification  et  de  1 ' analyse  de  1 ' erreur  systematique  et  des 
erreurs/incertitudes  de  modelisation  dans  1 1  analyse  technique. 

Les  sources  d ' erreur  systematique  et  d 1 erreurs/incertitudes  de 
modelisation  dans  1 ' analyse  de  fiabilite  basee  sur  les  MSEF  sont 
determinees  et  des  suggestions  sont  faites  en  vue  de 
1 1  application  de  formules  realistes  permettant  de  tenir  compte  de 
cette  classe  d 1  incertitudes .  Enfin,  on  fait  des  recommandations 
au  sujet  des  exigences  en  matiere  de  travaux  futurs  de 
developpement  qui  elargiront  la  portee  de  1 ' analyse  du  systeme 
STOVAST  et  permettront  des  applications  directes  aux  structures 
des  navi res . 
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1.  BACKGROUND 

1.1  General 

Issues  of  structural  safety  and  reliability  are  of  concern  to  design  engineers  as  well  as  those 
entrusted  with  the  responsibility  for  developing  design  criteria  and  operating/maintenance  schedules  and 
practices.  While  an  engineer’s  primary  task  in  designing  a  component  or  system  is  to  adequately  assure 
its  satisfactory  performance  during  its  useful  life,  the  planning,  design,  and  operation  of  most  (if  not 
all)'  engineering  systems  are  carried  out  under  conditions  of  uncertainty.  This  uncertainty  arises  due 
to  a  myriad  of  factors.  These  include  loads,  structural  strength  determined  by  such  uncertain  parameters 
as  geometry,  material  properties  and  sectional  properties,  boundary  conditions,  conditions  of 
manufacture,  fabrication,  construction,  storage,  and  environmental  conditions.  There  are  also 
uncertainties  associated  with  mathematical  modelling  which  usually  employ  some  simplifying 
assumptions,  the  limited  accuracy  of  the  numerical  methods  and  computing  machines  employed,  and 
unintended  human  errors.  A  comprehensive  description  of  various  types  of  uncertainty  that  may  be 
associated  with  structural  reliability  can  be  found  in  the  monograph  by  Melchers  [1]. 

Even  though  the  presence  of  these  uncertainties  is  well  recognized,  traditional  approaches 
to  structural  design  and  performance  assessment  have  treated  the  uncertain  parameters  as  deterministic 
constants  and  accounting  for  the  possible  variabilities  in  these  quantities  by  means  of  empirical  safety 
factors.  These  safety  factors  are,  in  many  cases,  wrongly  assumed  or  wholly  unknown  and  are  not 
founded  on  any  rational  philosophy.  Associated  with  a  given  nominal  safety  factor,  no  matter  how 
conservative  it  may  be,  there  is  invariably  some  underlying  probability  of  failure.  This  is  why 
engineering  structures  designed  with  the  safety  factor  approach  are  still  known  to  fail.  The  message 
from  this  is  that,  because  of  unavoidable  uncertainties,  the  assurance  of  structural  safety  or  performance 
cannot  be  absolute  and  so  there  is  always  invariably  some  underlying  probability  of  failure. 
Unfortunately,  conventional  deterministic  methods  of  analysis,  hereinafter  referred  to  as  the  traditional 
approach,  is  not  capable  of  quantifying  this  probability  of  failure. 

The  probabilistic  reliability  approach  is  proving  to  be  a  very  effective  and  rational 
procedure  for  accounting  for  the  uncertainties  inherent  in  any  engineering  system.  The  approach  has 
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1.2 


been  so  successful  in  the  characterization  of  uncertainties  and  the  prediction  of  structural  behaviour  in 
the  aerospace  and  nuclear  industries  that  it  has  formed  an  important  part  of  design,  analysis, 
maintenance,  and  inspection  planning  strategies.  Currently,  reliability-based  approaches  are  being 
developed  and  applied  by  NASA  in  conjunction  with  Southwest  Research  Institute  (SwRI)  to  the  design 
of  the  next  generation  Space  Shuttle  main  engine.  There  is  no  doubt  that  the  reliability  approach  is 
gaining  wider  acceptance  in  the  engineering  community,  especially  as  it  is  capable  of  assisting  in 
optimization  in  design  and  economization  in  maintenance  practices. 

1.2  Reliability  of  Shin  Structures 

In  the  same  manner  as  for  aerospace,  nuclear,  and  other  structures,  reliability  concepts  offer 
the  potential  for  achieving  greater  understanding  of  the  safety  and  performance  assessment  of  ship 
structures.  Although  the  most  obvious  source  of  uncertainties  in  the  case  of  ship  structures  is  the 
random  loading  induced  by  the  random  ocean  wave  environment,  there  is  no  doubt  that  the  structural 
systems  themselves  contain  inherent  uncertainties  in  material  properties,  section  properties,  structural 
geometry,  damping,  and  boundary  conditions.  Random  vibration  analysis  theory  is  applicable  to  cases 
in  which  only  the  loading  is  considered  to  be  random.  However,  when  the  variabilities  in  structural 
and  material  properties  are  significant  enough,  a  foil  scale  probabilistic  analysis  is  called  for. 

Mansour  [2]  was  the  first  to  put  forward  the  idea  of  probabilistic  design  concepts  in  marine 
structures  in  his  analysis  of  a  Mariner  ship.  Since  then,  the  field  has  expanded  enormously  with  the 
majority  of  applications  directed  to  offshore  structures.  A  review  of  reliability  methods  as  applied  to 
ship  structures  up  to  1980  was  given  by  Stiansen  et  al.  [3].  In  those  early  applications,  probabilistic 
analysis  and  reliability  methods  were  mainly  utilized  to  check  and  compare  the  safety  level  of  existing 
designs.  In  more  recent  times,  however,  some  civil  and  offshore  engineering  regulatory  bodies  and 
classification  societies  have  implemented  reliability  concepts  for  determining  their  code  requirements. 
They  are  also  increasingly  relying  on  the  use  of  probabilistic  methods  for  assessing  the  uncertainties 
underlying  the  design  or  performance  variables  and  for  selecting  the  associated  partial  safety  factors 
suitable  for  use  in  standards  and  codes.  A  comparison  of  various  approaches  to  estimating  reliability 
(ranging  from  first  order  methods  to  exact  (Level  III)  formulations)  was  made  by  Mansour  [4]  and 
applied  to  18  ocean  going  vessels. 
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Interest  in  ship  structural  reliability  is  increasing  as  evidenced  by  the  level  of  research 
activity  in  the  field.  The  American  Bureau  of  Shipping  (ABS)  and  the  Ship  Structure  Committee  (SSC) 
of  the  United  States  Coast  Guard  are  good  examples  of  agencies  that  are  continuously  working  on  a 
program  of  action  to  effect  applications  of  probabilistic  techniques  for  the  assessment  of  ship  structural 
reliability. 


Recently,  there  has  been  some  interest  in  developing  tools  for  reliability  assessment  of  ship 
structures  within  a  framework  that  permits  detailed  three-dimensional  analysis.  Just  like  its  deterministic 
counterpart,  the  stochastic  finite  element  method  is  the  most  viable  vehicle  for  achieving  this  objective. 
This  was  exemplified  in  a  recent  work  by  Moore  [5]. 

The  Defence  Research  Establishment  Atlantic  (DREA)  has  embarked  on  an  impressive 
research  and  development  program  on  the  applications  of  probabilistic  mechanics  technology  to  ship 
structures  in  the  last  four  years.  These  efforts  include  characterization  of  the  ocean  wave  environment, 
the  probabilistic  characterization  of  wave-induced  loads  within  the  framework  of  the  finite  element 
method  (FEM),  three-dimensional  finite  element  random  response  computations,  stochastic  fatigue 
damage  estimation,  and  applications  of  structural  reliability  theory. 

1.3  The  Stochastic  Finite  Element  Method 

Over  the  years,  the  fields  of  probabilistic  mechanics  and  computational  mechanics  have 
proceeded  vigorously  to  what  may  be  safely  referred  to  as  a  mature  stage  of  development.  Progress 
in  the  two  fields,  however,  has  proceeded  rather  independently  for  two  major  reasons.  First,  until 
recently,  and  except  for  a  few  exceptions,  research  in  each  of  the  above  fields  have  been  conducted 
by  separate  groups  that  are  essentially  mutually  exclusive.  Secondly,  the  application  of  stochastic 
methods  to  engineering  structures  generally  requires  considerable  computational  resources  and  resulted 
in  concerns  regarding  feasibility  of  the  tool  for  complex  models.  The  advent  of  powerful  computer 
hardware  and  the  development  of  novel  numerical  techniques  are  making  it  possible  to  address 
problems  on  a  more  realistic  basis.  Furthermore,  there  is  a  better  exchange  among  these  groups  of 
researchers  through  evaluation  and  verification  of  ideas  on  the  application  of  computational  methods 
to  problems  in  stochastic  mechanics. 
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The  field  of  computational  stochastic  mechanics  has  emerged  to  provide  the  marriage  of 
advanced  techniques  in  probabilistic  mechanics  and  computational  mechanics.  The  finite  element 
method  in  conventional  deterministic  analysis  has  gained  overwhelming  acceptance  in  the  mechanics 
community  because  of  its  robustness  and  suitability  for  the  modelling  of  complicated  structures.  For 
the  same  reason,  the  stochastic  finite  element  method  (SFEM)  is  emerging  as  a  very  powerful  tool  in 
probabilistic  structural  analysis.  The  basic  objective  of  SFEM  is  to  compute  the  probabilistic  structural 
response  in  terms  of  either  the  response  statistics  or  the  probability  of  failure  corresponding  to  a 
particular  performance  criterion  (or  set  of  criteria)  from  the  probabilistic  description  of  the  structure. 
The  latter  objective  is  especially  significant  from  a  practical  point  of  view  and  constitutes  the  basis  of 
stochastic  finite  element  based  structural  reliability  analysis.  Several  important  technical  issues  that 
are  of  direct  relevance  to  the  application  of  stochastic  finite  elements  for  structural  reliability  analysis 
were  discussed  in  a  recent  technical  note  by  the  first  author  [6]. 

1.4  Scone  of  the  Present  Work 

The  objectives  of  the  current  contract  work  consist  of  the  following: 

i.  -  A  state-of-the-art  literature  review  of  methodologies  for  stochastic  finite  element-based 

structural  reliability; 

ii.  Formulation  and  implementation  of  a  model  computer  program  for  demonstrating  the  best 
approaches  for  SFEM-based  reliability  analysis;  and 

iii.  Study  of  bias  and  modelling  errors  associated  with  SFEM-based  reliability  analysis. 

These  objectives  were  met  in  the  work  reported  in  the  present  document.  The  report  has 
seven  chapters.  Chapter  2  summarizes  the  literature  review  effort.  As  an  extension  of  the  literature 
review,  case  studies  of  two  existing  SFEM  based  computer  programs  are  given  in  Chapter  3.  The 
design  and  implementation  of  the  model  computer  program  (STOVAST)  for  performing  SFEM-based 
reliability  analysis  is  described  in  Chapter  4.  This  also  includes  the  development  of  the  random  field 
discretization  module  (RANFLD)  which  acts  as  a  pre-processor  to  program  STOVAST. 
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Some  demonstration/verification  problems  for  the  modei  code  are  presented  in  Chapter  5. 
In  Chapter  6,  a  discussion  on  bias/modelling  errors  associated  with  SFEM-based  reliability  analysis 
is  given.  Chapter  7  concludes  the  report  with  discussions  on  the  level  of  success  achieved  in  the 
project  and  recommendations  for  future  work  that  will  ensure  that  the  laudable  objectives  of  the  project 
will  be  vigorously  pursued. 
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2.  A  STATE-OF-THE-ART  REVIEW  OF  STOCHASTIC  FINITE  ELEMENT  BASED 

RELIABILITY  ANALYSIS  METHODOLOGIES _ 

2.1  Introduction 

In  this  chapter,  a  comprehensive  review  of  methodologies  for  stochastic  finite  element 
(SFEM)  based  reliability  analysis  is  presented.  This  review  covers  the  pertinent  subject  areas  which 
include  random  field  discretization,  general  approaches  for  SFEM-based  reliability  analysis,  reliability 
analysis  algorithms,  and  utilization  of  existing  FEA  computer  programs. 

2.2  Representation  and  Discretization  of  Random  Fields 

2.2.1  Definition 

A  random  field  is  a  family  of  random  variables  depending  on  more  than  one  deterministic 
parameter  and,  as  such,  it  may  be  viewed  as  a  multi-dimensional  random  process.  They  are  used  to 
represent  physical  systems  that  have  attributes  which  exhibit  complex  patterns  of  spatial  and  temporal 
variations.  Examples  of  random  fields  include  random  pressure  loading  due  to  ocean  waves,  stresses 
induced  by  turbulent  boundary-layer  pressure,  material  property  (such  as  Young’s  modulus)  with 
random  spatial  variations,  and  temperature  fields  in  the  main  engine  of  a  space  vehicle. 

Random  field  models  of  complex  stochastic  systems  serve  multiple  purposes.  They  provide 
a  format  for  efficient  description  of  the  random  variation.  They  provide  the  basis  for  predicting  system 
response  and  performance.  They  also  permit  the  assessment  of  the  impact  of  alternative  strategies  in 
diverse  situations  (such  as  in  system  design  or  data  requisition).  In  general,  it  is  usually  impractical  and 
uneconomical  to  sample,  for  example,  material  properties  of  a  random  medium  at  all  locations. 
Therefore,  the  tasks  of  prediction,  analysis,  and  decision-making  must  usually  proceed  on  the  basis  of 
incomplete  information  about  the  medium.  An  important  purpose  of  random  field  description  is  to 
permit  a  meaningful  representation  that  can  be  used  for  modelling.  A  detailed  exposition  of  the  theory 
of  random  fields  and  their  engineering  applications  is  given  in  the  monograph  by  Vanmarcke  [7], 
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2.2.2  Random  Field  Element  Discretization  Techniques 

2.2.2. 1  Objective  of  the  Discretization  Process 

The  terminology  "stochastic  finite  element  analysis"  (SFEA)  is  loosely  employed  to  describe 
the  explicit  stochastic  treatment  of  uncertainty  in  any  system  quantity  within  the  framework  of  the  finite 
element  method.  In  the  strict  sense,  however,  the  distinguishing  feature  of  stochastic  finite  element 
analysis  is  that  it  involves  the  discretization  of  the  parameter  space  of  a  random  field  of  material 
properties  and/or  loads. 

For  SFEM-based  structural  reliability  analysis,  in  particular,  discretization  is  required  to 
transform  the  spatial  distribution  of  the  field  to  "point"  values,  i.e.  random  variables.  The  resulting 
reliability  model  then  becomes  of  the  (generally  correlated)  random  variable  type,  for  which  the  vast 
majority  of  existing  computational  algorithms  for  reliability  analysis  were  designed  for.  In  this 
connection,  the  primary  statistics  of  the  random  variables  required  for  the  reliability  analysis  are  the 
mean  values  and  the  covariance  matrix  which  are  obtained  from  the  discretization  process. 

The  available  techniques  may  be  broadly  classified  into  three  groups.  In  one  group  of 
methods,  the  domain  of  the  field  is  discretized  into  a  mesh  of  random  field  elements  (not  necessarily 
coinciding  with  the  finite  element  mesh)  and  the  value  for  each  element  is  described  by  a  single  random 
variable.  The  spatial  averaging  methods,  the  midpoint  method,  the  nodal  point  method,  and  the 
interpolation  method  fall  into  this  category.  In  another  approach,  the  discretization  of  the  domain  of 
the  field  is  avoided,  and  instead  series  expansion  methods  are  used  to  model  the  field  as  a  series  of 
shape  functions  with  random  coefficients.  Finally,  the  third  group  of  methods  is  based  on  simulation. 
The  procedures  pertaining  to  each  of  these  groups  of  approaches  are  briefly  highlighted  in  what 
follows. 
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2.2.3  The  Spatial  Averaging  Method 

The  theory  of  local  averages  of  homogeneous  random  scalar  fields  over  rectangular  domains 
is  well  established.  A  detailed  discussion  on  the  subject  was  provided  by  Vanmarcke  in  his  textbook 
[7],  and  forms  the  basis  of  the  spatial  averaging  method. 

The  spatial  averaging  method,  suggested  by  Vanmarcke  and  Grigoriu  [8],  uses  the  local 
average  of  the  field  over  a  random  field  element  to  represent  the  random  quantity  for  the  element.  For 
a  random  field  W (X),  the  discretized  value  for  an  element  is  given  by: 

W.  =  —  /  W(X)  dQ  (2.1) 

Qi 

where  is  the  domain  of  the  element.  For  homogeneous  fields  and  rectangular  elements  whose  edges 
are  parallel  to  the  coordinate  axes,  expressions  for  the  covariances  of  the  discretized  variables  W;  in 
terms  of  the  autocovariance  function  of  W(X)  were  derived  by  Vanmarcke  [7], 

According  to  Der  Kiureghian  and  Ke  [9],  the  spatial  averaging  method  yields  accurate 
results  (even  for  coarse  meshes)  for  Gaussian  random  fields.  However,  it  has  two  major  shortcomings 
in  the  context  of  application  to  finite  element  reliability  analysis.  First,  for  a  two-dimensional  or  three- 
dimensional  continuum  of  arbitrary  shape  it  is  not  always  possible  to  discretize  the  domain  into 
rectangular  elements.  Several  approximation  schemes  that  have  been  suggested  in  the  literature  to  deal 
with  cases  involving  non-rectangular  elements  are  known  to  introduce  errors  in  the  computed 
covariance  matrix.  These  errors  may  lead  to  very  inaccurate  reliability  results  and  are  indeed  capable 
of  causing  the  FORM  and  SORM  algorithms  to  breakdown  in  some  cases  [10]. 

The  efforts  of  Zhu  and  coworkers  are  among  the  most  notable  attempts  to  relax  the 
restrictions  on  the  use  of  local  averages.  Zhu  and  Ren  [11]  extended  the  theory  to  deal  with 
homogeneous  and  rectangular  isotropic  random  vector  fields.  A  generalization  of  the  work  by  Zhu  and 
Wu  [12]  then  followed  in  which  the  condition  of  isotropy  was  relaxed  to  one  of  quadrant  symmetry. 
In  addition,  the  work  provided  the  capability  for  modelling  nonhomogeneous  random  vector  fields  and 
permitted  application  to  nonrectangular  domains.  The  allowance  for  irregular  domains  was  made 
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possible  by  the  use  of  Gaussian  quadrature  for  evaluating  the  means  and  covariances  of  local  averages. 
The  procedure  was  demonstrated  in  a  recent  work  involving  applications  to  SFEM  analysis  of  static, 
eigenvalue,  and  stress-intensity  factor  problems  [13]. 

Another  drawback  of  the  spatial  averaging  method  is  that  the  probability  distribution  of  W; 
is  difficult  or  impossible  to  obtain,  except  when  the  underlying  random  field  is  Gaussian,  in  which  case 
W;  is  also  Gaussian.  Essentially,  therefore,  the  use  of  this  method  is  restricted  to  Gaussian  fields. 


2.2.4  The  Midpoint  Method  and  the  Nodal-point  Method 


The  midpoint  and  nodal-point  methods  are  two  point  discretization  methods  which  represent 
the  uncertainties  of  a  random  field  by  the  values  at  some  specific  points.  In  the  midpoint  method  [9, 
14],  the  element  random  variable  is  defined  as  the  value  of  the  random  field  at  the  centroid  of  the 
element.  Hence,  the  randomness  in  a  random  field  element  i  is  represented  by  the  random  variable: 

Wj  =  W(X})  (2-2> 


in  which, 


NNODE 


NNODE  j.i 


2  X, 


(2.3) 


are  the  coordinates  of  the  centroid,  where  NNODE  is  the  number  of  nodes  of  the  random  field  element, 
and  X  *  are  the  nodal  coordinates. 


The  nodal-point  method  represents  the  random  field  in  terms  of  the  values  at  the  nodal 
points  of  the  finite  element  mesh.  In  this  method,  the  randomness  of  the  field  at  node  i  is  represented 
by: 

Wj  =  W(X.)  (2-4) 

where: 

Xj  «  Xd  (2-5> 


are  the  coordinates  of  node  i. 
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In  both  the  midpoint  method  and  the  nodal-point  method,  the  mean,  variance,  and  marginal 

distribution  of  W;  are  the  same  as  those  of  the  process  at  point  X{.  The  correlation  coefficient  matrix 
of  W;  is  directly  computed  in  terms  of  the  auto-correlation  coefficient  function  of  the  random  field, 

Pwiwj  =  Pww  CVj)  (2-6) 

and  the  joint  distribution  for  any  set  of  W;  is  given  by  the  specified  distribution  of  the  random  field. 

As  pointed  out  by  Der  Kiureghian  and  Ke  [9],  the  midpoint  method  tends  to  overrepresent 
the  variability  of  the  field  within  each  element,  and  it  does  not  provide  as  accurate  a  result  as  the 
spatial  averaging  method  for  a  coarse  random  field  mesh.  However,  these  point  discretization  methods 
have  three  advantages.  First,  no  complicated  computations  are  required  for  the  covariance  matrix  and 
the  method  is  easy  to  implement.  Second,  the  correlation  coefficient  matrix  obtained  by  Equation  (2.6) 
is  always  positive-definite,  provided  a  valid  autocorrelation  function  is  specified.  Hence,  the  numerical 
stability  problem  arising  in  the  spatial  averaging  method  does  not  exist  in  this  case.  Most  importantly, 
the  distribution  information  on  the  discretized  variables  W;  is  retained  and  the  method  is  not  restricted 
to  Gaussian  random  fields. 


2.2.5 


The  approach  employed  in  the  interpolation  method  suggested  by  Liu  et  al.  [15]  is  to 
represent  the  random  field  in  terms  of  an  interpolation  rule  involving  a  set  of  deterministic  shape 
functions  and  the  random  nodal  values  of  the  field.  Thus,  the  random  field  W (X)  is  discretized  into 
q  random  variables,  W;,  i  =  1,  ...,  q.  The  value  at  an  arbitrary  point  is  obtained  by  the  following 
interpolation  rule: 

W(X)  =  S  Ni(X)Wi  C2-7) 

i-l 

where  W-,  is  the  value  of  W(X)  at  node  X;,  and  N,(X)  are  shape  functions.  The  number  q  is  not 
necessarily  equal  to  the  number  of  finite  elements  and  the  shape  functions  N;(X)  need  not  be  the  same 
as  the  finite  element  interpolation  functions  for  the  displacement  field.  Since  the  choice  of  the  q  nodal 
points  and  shape  functions  is  arbitrary,  the  interpolation  method  constitutes  a  class  of  random  field 
discretization  methods.  In  particular,  it  is  easy  to  observe  that  if  the  nodes  are  chosen  to  be  the 
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centroids  of  the  random  field  elements  and  the  shape  functions  are  assumed  to  be  unity  inside  each 
element  and  zero  elsewhere,  the  interpolation  method  becomes  identical  to  the  midpoint  method 
described  in  Section  2.2.4. 

Liu  et  al.  [15]  further  suggested  a  method  to  reduce  the  number  of  random  variables  W-r 
The  random  vector  W  is  transformed  into  an  uncorrelated  random  vector  C  by: 

C  =  l|;TW  (2-8) 

such  that  the  covariance  matrix  of  C,  i.e.  Cov  (C,  C),  is  diagonal.  The  orthogonalization  matrix  \j/  is 
obtained  by  solving  the  eigenproblem: 

Cov(W,W)iJ;  =  i[rA  (2-9) 

where  A  is  the  eigenvalue  matrix  containing  the  variances  of  C.  Liu  et  al.  [15]  noted  that  a  good 
approximation  (within  the  context  of  second-moment  analysis)  of  the  random  field  can  be  obtained  by 
retaining  only  the  Q  with  large  variances,  thus  reducing  the  number  of  random  variables.  It  should  be 
emphasized,  however,  that  this  reduction  is  only  applicable  to  Gaussian  random  fields.  This  is  because 
the  distribution  of  C  is  generally  unknown  or  difficult  to  obtain  unless  W(X)  is  Gaussian. 

2.2.6  Series  Expansion  Methods 

Two  series  expansion  methods  have  been  suggested  for  second-moment  stochastic  finite 
element  analysis.  One  is  the  basis  random  variable  method  proposed  by  Lawrence  [16],  and  the  other 
is  the  kernel  expansion  method  proposed  by  Spanos  and  Ghanem  [17]. 

In  the  first  method,  the  random  field  is  expanded  into  a  double  series  in  the  form: 

W(X)  =  E  2  W..e.<j).(X)  (2-10) 

i=0  j-1  1  J 

in  which  j(X)  are  a  set  of  linearly  independent  shape  functions,  e;  are  independent  basis  random 
variables  having  the  properties: 
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E[ej] 


E[eiej]  =  6.. 


(2.11)a 
(2. 1  l)b 


in  which  5y  is  the  Kronecker  delta,  and  the  coefficients  Wg  are  determined  by  least-square  fitting  to  the 
moment  functions  of  the  function  field.  In  reference  [16],  Legendre  polynomials  were  suggested  for 
the  shape  functions,  0,(X).  In  applications,  only  a  few  dominant  terms  are  included  in  the  expansion. 

The  kernel  expansion  method  [17]  employs  the  Karhunen-Loeve  orthogonal  expansion  to 
decompose  a  one-dimensional  random  field.  The  random  field  is  expanded  into  the  sum  of  its  mean 
function  and  a  single  series: 


W(X)  =  pw(X)  +  2  Wj  4>i(X) 

i-0 


(2.12) 


where  W{  are  random  coefficients  independent  of  X,  and  \  and  ^(X)  are  the  eigenvalues  and 
eigenfimctions  of  the  covariance  kernel,  respectively.  The  latter  are  obtained  as  the  solutions  of  the 
eigenvalue  problem: 


I  COV(X,t)(|)i(X)dt  =  AjCjjjC X) 


(2.13) 


The  series  in  Equation  (2.12)  are  truncated  after  the  first  few  dominant  terms  just  as  in  the 
basis  random  variable  method.  Since  this  series  has  a  zero  mean  and  the  eigenfunctions  are  orthogonal, 
the  random  coefficients  W;  have  properties  similar  to  the  basis  random  variables,  i.e.  Equations  (2. 1  l)a 
and  (2. 1  l)b.  One  major  obstacle  of  this  method  is  the  difficulty  in  solving  the  eigenvalue  problem  in 
Equation  (2. 13)  for  arbitrary  geometry  and  boundary  conditions.  As  noted  by  Liu  and  Der  Kiureghian 
[10],  this  would  be  a  particularly  difficult  task  if  the  method  is  to  be  extended  to  two-dimensional  or 
three-dimensional  random  fields.  Details  of  this  method  are  extensively  discussed  and  exemplified  in 
a  recent  monograph  by  Ghanem  and  Spanos  [18]. 


Both  series  expansion  methods  described  above  are  strictly  applicable  only  to  Gaussian 
random  fields  by  virtue  of  the  central-limit  theorem.  Therefore,  they  are  appropriate  for  second- 
moment  analysis,  or  for  reliability  analysis  when  the  random  fields  are  truly  Gaussian. 
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2.2.7  Simulation  Methods 

Digital  simulation  is  a  technique  of  last  resort  to  solve  complex  engineering  problems 
involving  random  variation.  Assuming  a  deterministic  solution  method  is  available,  the  simulation 
technique  consists  of  repeatedly  generating  sets  of  observations  of  all  the  random  variables  or  functions 
involved  in  the  calculation,  solving  the  deterministic  problem  associated  with  each  set  of  the 
observations,  and  evaluating  appropriate  statistics  of  response  and  performance  measures. 

Shinozuka  has  developed  and  advocated  the  use  of  Monte  Carlo  simulation  techniques  in 
the  field  of  engineering  mechanics  over  the  last  two  decades.  A  synopsis  of  his  works  in  this  regard 
is  well  documented  in  a  recent  review  article  co-authored  with  Deodatis  [19].  In  particular,  Yamazaki 
and  Shinozuka  [20]  proposed  an  interactive  procedure  to  simulate  non-Gaussian  stochastic  fields  and 
later  introduced  statistical  preconditioning  to  reduce  the  sample  size  in  another  work  [21]. 

Recently,  Fenton  and  Vanmarcke  [22]  presented  a  method  of  generating  realizations  of  a 
discrete  "local  average"  Gaussian  random  process  using  a  local  average  subdivision  technique.  One  of 
the  advantages  cited  by  the  authors  is  the  suitability  for  finite  element  models  that  employ  efficient  low- 
order  interpolation  functions  in  which  each  local  average  becomes  an  element  property.  No  structural 
applications  of  this  technique  has  been  reported  in  the  open  literature. 

More  recent  contributions  in  connection  with  the  application  of  simulation  schemes  for 
modelling  and  discretization  of  random  fields  include  the  works  of  Bielewicz  et  al.  [23]  on 
nonhomogeneous  scalar  fields,  Grigoriu  [24]  on  the  application  of  the  sampling  theorem  and  Poirion 
[25]  on  the  simulation  of  non-Gaussian  fields. 

A  major  disadvantage  of  simulation  schemes  is  the  enormous  amounts  of  computation  times 
associated  with  their  application.  For  a  large  scale  structure  like  a  ship,  therefore,  these  procedures  are 
not  recommended. 
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2.2.8  Other  Methods' for  Random  Field  Discretization 

A  new  technique  was  recently  introduced  for  random  field  discretization  by  Li  and  Der 
Kiureghian  [26].  The  technique  is  referred  to  as  the  optimal  linear  estimation  (OLE)  method.  For  a 
random  field  defined  in  a  domain  Q,  a  linear  function  of  the  nodal  values  W(X;)  is  proposed  in  the 
form: 

N 

W(x)  =  a(x)  +  Sbj(x)T  v(Xj)  ,  xeQ  (2-14) 

i-1 


where  N  denotes  the  number  of  nodal  points  in  the  domain.  The  functions  a(x)  and  b;(x)  are 
determined  optimally  by  minimizing  the  variance  of  the  error  W(x)  -  W(x).  The  minimization  is 

A 

subject  to  the  condition  that  W(x)  is  an  unbiased  estimator  in  the  mean,  that  is. 


E[W(x)  -  W(x)]  =  0  . 

This  results  in  the  representation: 

W(x)  =  p(x)  +  S^oow  ^ww  (W-p.),xeQ 


(2.15) 


(2.16) 


in  which  W  denotes  the  vector  of  nodal  values  (i.e.  the  random  variables)  and  £W(X)W  denotes  an  Nxl 
vector  containing  the  covariances  of  W(x)  with  the  elements  of  W. 


The  OLE  method  is  believed  to  be  always  superior  to  the  midpoint,  spatial  averaging,  and 
shape  function  methods,  and  its  efficiency  can  be  further  improved  by  using  eigenvalue  expansion. 
The  OLE  method  can  also  be  applied  for  non-Gaussian  random  fields  by  utilizing  the  Nataf  model  of 
Liu  and  Der  Kiureghian  [27].  The  Nataf  model  assumes  that  a  transformed  process  z(x)  can  always 
be  found  such  that 

W(x)  =  Q-1  p(z(x))]  ,  (2-17> 

where  4?  is  the  standard  normal  cumulative  distribution  function.  The  key  to  the  superior  accuracy  and 
efficiency  of  this  method  lies  in  the  utilization  of  shape  functions  that  take  the  correlation  structure  of 
the  random  field  into  account.  The  technique,  therefore,  appears  to  be  a  very  promising  tool  in 
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stochastic  finite  element  analysis  where  it  is  always  desirable  to  represent  a  random  field  with  as  few 
random  variables  as  possible. 

2.3  Methodologies  for  SFEM-Based  Structural  Reliability  Analysis 

2.3.1  Merit  of  SFEM-Based  Approach 

In  a  stochastic  finite  element  based  approach  to  structural  reliability  analysis,  the  methods 
of  structural  reliability  provide  the  basis  for  modelling  and  analysis  of  uncertainties  and  computation 
of  probabilities,  while  the  finite  element  method  provides  the  necessary  computational  framework  for 
analyzing  complex  structures.  The  application  of  these  two  concepts  in  a  closely  integrated  fashion  is 
the  basis  of  SFEM-based  reliability  analysis,  and  results  in  a  powerful  tool  for  realistically  dealing  with 
practical  engineering  problems. 

The  merit  of  the  SFEM-based  approach  was  clearly  demonstrated  by  Orisamolu  et  al.  [28] 
in  a  recent  Martec  study.  In  that  study,  one  of  the  most  powerful  and  internationally  recognized 
general  purpose  commercial  probabilistic  analysis  programs  (PROBAN®)  was  applied  for  the  systems 
reliability  analyses  of  the  two-storey  braced  frame,  and  then  the  tubular  steel  jacket  offshore  platform 
whose  finite  element  models  are  shown  in  Figures  2.1a  and  2.2a,  respectively.  The  external  loads  (as 
shown  in  Figures  2.1b  and  2.2b)  and  some  other  strength  and  modelling  parameters  were  considered 
to  be  stochastic.  Details  of  the  description  and  analysis  of  this  problem  are  given  in  reference  [28].  The 
offshore  platform  problem  is  discussed  further  below  since  it  sufficiently  illustrates  the  merit  of  the 
integrated  approach. 

For  the  offshore  platform  problem,  three  failure  modes  were  considered  for  the  two- 
dimensional  tubular  members.  These  include  yielding,  stability,  and  punching.  Performance  functions 
corresponding  to  these  failure  modes  are  described  below. 

The  performance  function,  gy,  corresponding  to  yielding  failure  of  a  two-dimensional 
tubular  element  is  given  by: 
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gy(X)  =  zt 


JlT  _  ML 
nfJ  mf 


(2.18) 


where  Z:  is  a  modelling  uncertainty  variable.  N  is  the  axial  force  on  the  element  defined  as: 


N  = 


NLOAD 

2 


j-1 


a.P. 

j  j 


(2.19) 


where  NLOAD  is  the  number  of  independent  loads  acting,  a:-  are  force  coefficients  of  influence  and  P„ 
P 2 ,...,  P Nload  are  the  external  loads  modelled  as  stochastic  variables.  Figure  2c  shows  the  section 
forces  used  in  modelling  this  type  of  failure  element.  The  two  end-moments  Mt  and  M2  are  determined 
from: 


NLOAD 
Mj  =  2 

j-i 


b.P. 
u  i 


,  i=l,2 


(2.20) 


where  by  are  moment  coefficients  of  influence.  The  resultant  moment  in  the  actual  end  is,  for  end  I, 


M  =  Mx 


(2.21) 


and  for  end  J, 

M  =  Mj  (2-22) 

NF  in  Equation  (2.18)  is  determined  from  the  relation: 

Np  =  Yu  (d2  -  (d-t)2)/4  (2-23) 

while  MF  is  given  by 

MF  =  Y  (d3  -  (d-t)3)/6  (2-24) 

where  Y  is  the  yielding  stress,  d  is  the  tubular  diameter  and  t  is  the  tubular  thickness.  This  failure  is 
assumed  to  have  a  brittle  behaviour,  i.e.  it  is  assumed  that  the  element  is  removed  by  failure. 


For  stability  failure,  the  performance  function,  g„  for  a  two-dimensional  tubular  element 


is  given  by: 
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where: 


s.(X)  -  z,  -  ^ 


M 


(2.25) 


NF  =  YA  =  Yudt 
MF  =  Yd  h 


fYA  (1-  A.2/4) 

,  X  <;  s/2 

[YAA2 

,  X  >  \J2 

AY 

X2 

(2.26) a 

(2.26) b 

(2.26) c 

(2.26) d 


The  buckling  parameter  (X)  is  determined  by  the  expression: 


(2.27) 


In  the  above  relations,  Y  is  the  yielding  stress,  2^  is  a  model  uncertainty  variable  for  this  mode  of 
failure,  Z3  is  a  model  uncertainty  variable  associated  with  modelling  effective  column  lengths,  A  is  the 
section  area,  d  is  the  tubular  diameter,  t  is  the  tubular  thickness,  Y  is  the  cross-section  radius  of 
gyration,  L  is  the  column  length,  and  E  is  the  Young’s  modulus.  Again,  the  section  forces  used  in 
modelling  this  failure  element  are  as  shown  in  Figure  2.2c.  The  axial  force  N  and  the  two  moments 
Mj  and  M2  are  determined  by: 
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N 

NLOAD 

=  2  Vi 

(2.28)a 

NLOAD 

j=i 

(2.28)b 

NLOAD 

=  2,  V) 

J  =  1 

(2.28)c 

where  by,  b2j,  and  b3j  are  coefficients  of  influence  and  P,,  P2,...,Pnload  are  the  external  loads  modelled 
as  stochastic  variables.  The  moment  M  is  defined  as: 

M  =  Max  (Mp  Mj)  (2.29) 


The  punching  failure  element  has  the  performance  function,  gp,  given  by: 


gP(X)  =  Z4  - 


N„ 


ZjNu 


M 


V.2 


(2.30) 


where: 


YT2 

N  =  1  (3.4  +  19  B)ix 

(2.3  l)a 

sinQ 

-  YT2  0.8d  (3.4  +  19B)u 

(2.3  l)b 

sine 

|N  | 

[i  -  min  (1,  1.22  -  cl  ) 

1  '  7  att  r\m/ 

(2.3  l)c 

2YttDT 


in  which  Y  is  the  yielding  stress,  and  Z4  ,  Z5,  and  Z6  are  modelling  uncertainty  variables.  This  failure 
element  models  punching  failure  of  the  plane  K-joint  shown  in  Figure  2.2d.  It  is  assumed  that  the 
punching  failure  of  each  branch  can  be  considered  separately.  The  brace  has  a  diameter  d  and  thickness 
t,  while  the  chord  has  a  diameter  D  and  thickness  T.  The  angle  between  the  brace  and  the  chord  is  0. 
These  parameters  defining  the  geometry  of  the  joint  are  shown  in  Figure  2. 2d.  The  axial  force  in  the 
brace,  Nb,  and  in  the  chord,  Nc,  are  given,  respectively,  as: 
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NLOAD 

Nb  = 

E 

i=i 

aibpj 

(2.32) 

NLOAD 

(2.33) 

Nc  - 

E 

aicpi 

The  moment  in  the  actual  branch  is  given  as: 

M  = 

NLOAD 

E 

blpj 

(2.34) 

i-i 

In  the  application  of  PROBAN  for  the  probabilistic  analysis  of  this  offshore  platform,  each 
of  the  above  failure  modes  had  to  be  coded  and  linked  into  the  PROBAN  limit  state  function  library. 
It  can  be  seen  from  the  above  limit  state  equations  that  the  force  and  moment  influence  coefficients  are 
involved.  For  this  complicated  structure,  these  coefficients  are  best  determined  by  performing  a  linear- 
elastic  finite  element  analysis.  For  this  case  the  Martec  FEA  code  (VAST)  was  used  to  determine  these 
influence  coefficients.  They  are  to  be  supplied  as  part  of  the  input  data  for  each  of  the  components  that 
makeup  the  structural  system  in  the  description  of  the  probabilistic  model  to  PROBAN.  Now,  for  this 
offshore  platform,  a  total  of  about  62  failure  elements  consisting  of  36  yielding  elements,  18  stability 
elements,  and  eight  punching  elements  are  involved.  If  all  these  62  failure  elements  are  to  be  accounted 
for  in  the  evaluation  of  systems  reliability,  there  is  a  lot  of  tedium  associated  with  the  process  of 
preparing  the  enormous  amount  of  input  data  required.  This  is  especially  so  since,  for  the  punching 
failure  function,  the  calculation  of  angles  between  the  members  (which  is  not  an  output  of  the  FEA 
results)  is  required.  All  these  efforts  are  required  because  PROBAN  does  not  include  finite  element 
analysis.  The  computer  program  is  aimed  at  the  analysis  of  problems  whose  probabilistic  models  are 
well  defined  in  terms  of  explicit  limit  state  functions  and  statistical  characteristics  of  the  random 
variables. 


A  fully  integrated  and  coupled  finite  element  reliability  analysis  program  is  the  most  suitable 
for  the  analysis  of  the  offshore  platform  problem  described  above,  and  even  more  so  for  continuum 
structures.  For  such  an  integrated  program,  the  statistical  properties  of  the  random  variables  are  merely 
supplied  as  additional  information  to  the  finite  element  input  data  and  the  probabilistic  analysis  proceeds 
automatically.  The  user  is,  therefore,  alleviated  of  the  tedium  involved  in  the  detailed  description  of 
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the  probabilistic  model.  In  addition,  the  risk  of  errors  in  defining  the  input  data  is  considerably 
reduced.  For  example,  the  calculation  of  the  influence  coefficients  will  be  unnecessary  as  this  operation 
will  be  internally  taken  care  of. 

Furthermore,  in  the  use  of  SFEM-based  reliability  analysis  programs  for  problems  involving 
semi-ductile  or  ductile  elements  (such  as  the  stability  failure  element  described  above)  the  redefinition 
of  the  loads  is  automatically  accounted  for.  The  use  of  semi-ductile  or  ductile  elements,  in  strictly 
probabilistic  analysis  programs  involves  the  tedious  operation  of  analysis,  input  data  definition,  then 
re-analysis  and  redefinition  of  input  data,  and  so  on.  Also,  for  continuum  structures  (or  even  for 
discrete  structures  subjected  to  distributed  random  loads)  the  calculation  of  the  random  loads  acting  at 
several  locations  is  best  achieved  via  deterministic  or  stochastic  finite  element  discretization  procedures. 
Also,  this  scheme  is  indispensable  if  the  tubular  members  were  to  have  stochastic  spatial  variability  in 
the  Young’s  modulus  or  yield  stress.  Thus,  it  is  obvious  that  problem  formulation  and  analysis  can  be 
more  efficiently  performed  within  the  framework  of  an  integrated  finite  element  reliability  analysis 
computer  package. 

The  merit  of  SFEM-based  reliability  analysis  is  well  recognized  in  the  research  community. 
As  such,  it  represents  a  very  active  area  of  investigation.  Currently  available  methodologies  may  be 
classified  into  four  major  groups.  These  include  the  perturbation-based  SFEM  approach,  the  reliability- 
based  SFEM  approach,  response  surface  methods,  and  Monte  Carlo  simulation  methods.  The  weighted 
integral  method  has  also  recently  been  introduced  [29,  30]  to  avoid  a  direct  discretization  of  the  random 
fields  during  a  stochastic  finite  element  analysis.  The  basic  features  pertaining  to  each  of  these  groups 
of  techniques  are  briefly  highlighted  in  what  follows.  With  the  exception  of  the  Monte  Carlo  simulation 
and  the  weighted  integral  techniques,  these  methods  are  generally  based  on  the  combination  of 
stochastic  finite  element  discretization  techniques  and  the  well  known  first  order  reliability  methods 
(FORM)  and  second  order  reliability  methods  (SORM). 

2.3.2  Perturbation-Based  SFEM  Approach 

The  perturbation  approach  to  probabilistic  structural  analysis  was  introduced  about  two 
decades  ago.  Initial  applications  were  directed  at  the  study  of  the  eigenvalue  problem  related  to  the  free 
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vibration  of  structures  with  stochastic  mass  and  stiffness  matrices  and  the  solutions  of  linear  static 
problems  involving  loading  and  system  stochasticity.  The  work  of  Hisada  and  Nakagiri  [14]  represents 
one  of  the  modern  applications  of  this  approach  to  structural  safety  and  reliability  analysis.  In  that 
work,  SFEM  was  applied  for  the  evaluation  of  the  reliability  index  and  design  point  within  the 
framework  of  the  Advanced  First-Order  Second  Moment  (AFOSM)  method. 

Consider  the  linear  finite  element  equation 

KU  =  F 

where  K  is  the  global  stiffness  matrix,  U  is  the  vector  of  nodal  displacement  and 
load  vector.  These  stochastic  quantities  can  be  expressed  as 

K  =  K0  +  AK 

U  =  U0  +  AU 

F  =  F0  +  AF 


(2.35) 

F  is  the  global  nodal 

(2.36) a 

(2.36) b 

(2.36) c 


where  K0,  U0  and  F0  are  considered  to  be  deterministic  parts  and  the  increments  AK,  AU  and  AF  are 
considered  to  be  the  stochastic  parts  of  K,  U,  and  F,  respectively.  Substituting  Equations  (2.36)a-c 
into  (2.35)  gives 

[K0  +  AK]  [U0  +  AU]  =  p0  +  AF]  .  (2-37) 

on  neglecting  the  product  (AKAU)  and  separating  the  deterministic  and  stochastic  parts  of  Equation 
(2.37)  gives: 


K0U0  -  F0 


(2.38) 


and 

K0AU  -  AF-AKU0  (2-39) 

Equation  (2.38)  gives  the  finite  element  solution  at  the  deterministic  expansion  point.  From  the 
solution  of  Equation  (2.39),  the  second-order  variation  of  the  response  may  be  computed. 


P152713.PDF  [Page:  39  of  174] 


2.17 


A  more  rigorous  and  general  formulation  of  the  perturbation  approach  can  be  constructed 
using  Taylor  series  expansion.  This  also  paves  the  way  for  higher-order  approximations.  The  stiffness 
matrix  in  Equation  (2.35)  may  be  expanded  about  a  deterministic  stake  as: 

K  =  K0  +  S  Kj'Cx-Xj)  +  i.  2  2  K '(x-XjXx-x.),  (2-40) 

i=l  2  i«l  j»l] 

where  K;'  and  Ky'  are  the  first-order  and  second-order  partial  derivatives  of  the  stiffness  matrix  with 
respect  to  the  basic  variables  X.  The  displacement  and  load  vectors  in  Equation  (2.35)  may  be 
expanded  in  a  similar  fashion  and  the  response  computed  as  previously  described. 

When  the  mean  state  is  chosen  as  the  expansion  point  the  approach  is  referred  to  as  the 
mean-centred  perturbation  approach.  Using  mean-centred  perturbation  results  in  the  computation  of 
reliability  indices  according  to  the  First-Order  Second  Moment  (FOSM)  approach  of  structural 
reliability  analysis.  In  the  AFOSM,  however,  the  performance  function  is  expanded  not  about  the  mean 
values  of  the  basic  variables,  but  about  the  most  probable  failure  point.  Hisada  and  Nakagiri  utilized 
this  approach  and  also  presented  a  second-order  perturbation  formulation  in  reference  [14].  A  notable 
feature  of  the  formulations  is  that  the  stiffness  matrix  is  inverted  only  once  in  contrast  to  simulation 
or'response  surface  methods  in  which  many  inversions  of  the  stiffness  matrix  are  required.  The  key 
to  successful  solution  using  the  perturbation  approach  is  the  ability  to  compute  and  assemble  partial 
derivative  matrices  for  stiffness,  displacements,  and  loads.  Second-order  approximations  are  obviously 
more  accurate  than  the  first-order  approximations;  however,  these  involve  the  computation  and 
assembly  of  second-order  partial  derivative  matrices. 

Numerous  applications  of  the  perturbation  approach  have  been  reported  in  the  open 
literature.  Prominent  in  this  connection  are  the  works  of  Liu  et  ai.  [15,  31,  32],  in  which  applications 
were  investigated  for  linear  and  nonlinear  structural  dynamics  and  a  variational  formulation  of 
probabilistic  finite  elements  established. 

Although  the  formulation  of  the  perturbation  approach  is  mathematically  elegant,  its 
application  to  reliability  analysis  has  several  disadvantages.  The  mean-centred  perturbation  method 
suffers  from  the  well  known  invariance  problem  associated  with  FOSM.  Furthermore,  the  perturbation 
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methods  do  not  use  the  distribution  information  about  the  basic  random  variables,  even  if  it  is 
available.  This  is  a  serious  limitation,  unless  for  the  exceptional  cases  in  which  all  the  variables  are 
normally  distributed.  The  method  is  also  not  capable  of  producing  accurate  results  when  there  are 
large  variations  in  the  random  variables  defining  a  problem. 

2.3.3  The  Reliabilitv-Based  SFEM  Approach 

Current  methods  for  reliability  analysis  compute  the  reliability  index  /3  by  solving  the  limit 
state  equation  g(X)=0  explicitly.  However,  search  algorithms  that  do  not  rely  on  an  explicit  solution 
of  the  limit  state  equation  are  available.  These  algorithms  only  need  the  value  and  the  gradient  of  the 
performance  function  at  each  iteration  point. 

The  reliability  approach  has  been  formulated  and  applied  for  several  structural  problems 
by  Der  Kiureghian  and  his  co-workers.  In  [9,  33]  the  first-order  reliability  method  was  used  for  static 
analysis  of  linear  structures  with  random  properties  and  in  [16]  for  dynamic  analysis.  In  a  more  recent 
work  [10],  a  general  framework  for  finite  element  reliability  analysis  based  on  FORM  and  SORM  was 
presented.  New  expressions  for  the  required  gradients  of  the  response  of  geometrically  nonlinear 
structures  were  derived  and  implemented.  This  work  represents  the  first  application  of  the  finite 
element  reliability  method  (FERM)  in  conjunction  with  SORM  with  non-Gaussian  random  fields,  and 
with  system  reliability  analysis.  Arnbjerg-Nielsen  and  Bjerager  [35]  and  Mahadevan  [36]  has  also 
developed  and  implemented  a  reliability-based  SFEM  approach  and  applied  it  to  the  modelling  of  frame 
structures. 


The  reliability  approach  has  a  significant  advantage  over  the  perturbation  approach  in  that 
information  about  the  distribution  of  the  random  variable  is  used.  Furthermore,  in  the  reliability 
approach,  the  probability  density  function  of  the  response  variable  (not  just  second  moment  statistics) 
can  be  obtained.  Since  the  computation  of  response  gradients  is  a  key  operation  in  the  implementation 
of  this  procedure,  the  use  of  the  adjoint  method  has  been  recommended  for  this  operation  [10,  37]. 
This  technique  is  estimated  to  be  capable  of  reducing  computation  times  by  a  factor  equal  to  the 
number  of  random  variables. 
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The  reliability-based  approach  to  SFEM  structural  reliability  analysis  is  especially 
compatible  with  the  algorithmic  structure  of  existing  FEA  codes  and  is  about  the  best  strategy  of  all 
presently  available  methodologies. 

2.3.4  Response  Surface  Methods 

The  response  surface  method  is  a  classical  statistical  technique  in  which  a  complex 
(computer)  model  is  approximated  by  a  simple  functional  relationship  between  the  output  quantities  and 
the  input  (basic)  variables.  The  approximation  is  usually  based  on  polynomial  functions  and,  often, 
linear  or  quadratic  response  functions  are  applied.  Adopting  the  simpler  response  functions  allows  an 
■  efficient  repeated  computation,  for  example,  as  may  be  needed  in  simulations  or  parameter  studies  in 
structural  reliability  analysis.  This  is  because  the  approximation  to  the  response  surface  rather  than  the 
original  limit  state  function  is  used  in  the  calculation  of  failure  probabilities. 

The  concept  of  response  surface  methods  has  been  used  when  approximating  costly  to 
compute  and/or  non-differentiable  limit  state  functions.  Within  the  framework  of  the  stochastic  finite 
element  method,  the  steps  required  for  the  implementation  of  the  response  surface  technique  were 
described  by  Favavelli  [38],  This  involves  the  application  of  regression  analysis  to  obtain  the 
polynomial  coefficients  involved  in  the  representation  of  the  limit  state  function  using  the  results  of 
several  numerical  experiments. 

The  explicit  representation  of  the  limit  state  function  g(x),  for  the  quadratic  approximation 
for  example,  takes  the  form: 

g(x)  =  a  +  SbjXj  +  2  2  C..X.X,  (2-41) 

i-1  i-1  j-1  J  J 

where  n  is  the  number  of  basic  random  variables  (xj  and  the  coefficients  a,  b,  and  c  are  to  be 
determined  from  numerical  experiments. 

The  works  of  Schueller  et  al.  [39],  Bucher  and  Bourgund  400  and  Bohm  and  Bruckner-Foit 
[41]  are  among  recent  efforts  at  promoting  the  application  of  response  surface  methods.  The  work  of 
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Bohm  and  Briickner-Foit  [41],  in  particular,  introduced  a  special  lack  of  fit  measure  and  formulated 
criteria  for  accepting  response  surface  models  in  structural  reliability  analysis.  Ghanem  and  Spanos 
[42]  proposed  a  Galerkin-based  response  surface  approach  in  which  the  surface  is  approximated  by  its 
projection  onto  a  complete  set  of  polynomials  that  are  orthogonal  to  the  Gaussian  measure.  These 
polynomials  are  known  as  polynomial  chaos  functions  and  are  believed  to  be  capable  of  yielding 
accurate  approximations  of  the  response  surface. 

So  far,  no  general  scheme  has  been  developed  to  efficiently  establish  linear  and  quadratic 
response  surfaces  for  reliability  computations.  Further  research  is  needed  to  establish  general,  efficient, 
and  robust  response  surface  methods  for  reliability  analysis.  Nevertheless,  this  methodology  appears 
to  be  a  promising  tool  for  large  scale  structures  as  recently  demonstrated  in  a  study  dealing  with  the 
fatigue  reliability  of  a  container  ship  structure  [43]. 

2.3.5  Monte  Carlo  Simulation  Methods 


The  direct  Monte  Carlo  simulation  method  was  used  in  many  early  works  in  stochastic  finite 
element  analysis.  In  this  brute  force  method,  deterministic  analysis  is  carried  out  for  a  series  of 
parameters  generated  in  accordance  with  their  probability  distribution.  The  desired  statistics  of  the 
response  quantities,  such  as  the  mean,  variance,  and  exceedance  probabilities,  are  then  evaluated  based 
on  the  generated  sample.  Applications  of  this  procedure  can  be  found  in  Vanmarcke  et  al.  [44]  and 
Takada  [45]. 

The  Monte  Carlo  simulation  method  has  the  advantage  that  it  is  adaptable  to  all  types  of 
problems  and  the  results  can  be  obtained  to  desired  accuracy.  However,  for  practical  problems  with 
many  random  variables  or  small  failure  probabilities  this  procedure  is  usually  too  expensive,  since  a 
large  number  of  solutions  are  needed  to  obtain  reliable  results.  Shinozuka  and  his  co-workers  [46] 
have  introduced  the  Neumann  (Monte  Carlo  simulation)  expansion  technique.  Computation  time  for 
this  technique  is  reduced  significantly  since  only  the  mean  stiffness  matrix  needs  to  be  decomposed  with 
this  formulation.  Other  schemes  have  also  been  proposed  to  improve  the  efficiency  of  the  simulation 
method.  However,  for  ship  structures,  this  procedure  is  not  recommended. 
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2.3.6  The  Weighted  Integral  Method 

Under  certain  conditions,  discretization  of  the  random  fields  in  SFEM  may  be  avoided 
through  the  application  of  the  weighted  integral  method.  This  method  was  recently  introduced  by 
Deodatis  [29].  It  was  applied  by  Deodatis  and  Shinozuka  [30]  for  the  calculation  of  the  response 
variability  and  reliability  of  stochastic  frame  structures  and  the  calculation  of  the  response  variability 
of  two-dimensional  stochastic  systems  by  Deodatis  et  al.  [47], 

Basically,  the  method  assumes  that  the  elastic  constitutive  matrix  (D)  in  the  element  stiffness 
formulation  can  be  expressed  as  a  product  of  a  scalar  random  field  E(x)  and  a  deterministic  matrix  D„, 
that  is: 

D  =  E(x)Do  (2.42) 

Then,  due  to  the  polynomial  nature  of  the  strain  matrix,  B,  the  element  stiffness  matrix  can  be  written 

as: 

K  =  S  X.K  ,  (2-43) 

e  i  a 

i-1 

where  n  denotes  the  number  of  distinct  polynomial  terms  in  the  coordinates  x  in  the  matrix  product 
BTB,  Kei  are  deterministic  matrices,  and  X;  are  random  variables  defined  by: 

X.  =  / P,(x)  E(x)  dQe  (2.44) 

0. 

in  which  P;(x)  denotes  the  i-th  distinct  polynomial  term  in  the  coordinates  x.  The  random  variables  X; 
are  interpreted  as  "weighted  integrals"  for  all  elements  and  completely  define  the  SFEM  problem 
without  the  need  to  discretize  the  random  fields. 

The  major  advantage  of  this  procedure  has  been  cited  as  the  circumvention  of  the  restriction 
(in  random  field  discretization)  that  the  finite  element  size  has  to  be  a  fraction  of  the  correlation 
distance  of  the  stochastic  field  involved  in  the  problem.  This  restriction  makes  necessary  the  use  of  a 
fine  mesh  to  accurately  describe  stochastic  fields  characterized  by  short  correlation  distances. 
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The  method,  however,  has  two  main  disadvantages.  First,  the  assumption  of  a  product  form 
assumed  in  Equation  (2.42)  is  generally  difficult  to  realize  in  practice.  Secondly,  the  method  assumes 
a  Gaussian  random  field  which  is  a  serious  drawback  in  reliability  analysis  where  probabilities  are 
sensitive  to  tails  of  the  underlying  distributions. 

2.3.7  Utilization  of  Commercial  FEA  Programs 

A  few  applications  involving  the  use  of  existing  general  purpose  commercial  finite  element 
programs  for  reliability  analysis  have  been  reported.  Gopalakrishna  and  Donaldson  [48]  utilized 
ANSYS  within  the  framework  of  the  advanced  first-order  second  moment  (AFOSM)  method  to  estimate 
failure  probabilities.  The  work  took  advantage  of  the  response  surface  based  design  optimization 
routine  already  available  in  ANSYS.  The  solution  of  the  reliability  problem  is  constructed  such  that 
the  basic  random  variables  represent  the  design  variables,  the  limit  state  function  constitutes  the 
constraint  function  and  the  reliability  index  is  treated  as  the  objective  function. 

Riha  et  al.  [49]  also  presented  a  coupling  of  a  fast  probability  integration  scheme  with  the 
finite  element  code  MSC/NASTRAN.  Again,  in  this  work,  advantage  was  taken  of  the  capability  for 
design  sensitivity  analysis  (DSA)  available  in  the  finite  element  system  for  the  computation  of  structural 
response  derivatives.  The  DSA  capability  however,  is  not  applicable  to  some  boundary  condition  and 
load  variables  and  so  cannot  be  used  for  the  calculation  of  the  partial  derivatives  of  such  random 
variables.  Provision  was  made,  however,  for  the  use  of  the  finite  difference  scheme. 

The  probabilistic  analysis  capabilities  described  in  the  preceding  couple  of  paragraphs  do 
not  fall  into  the  general  framework  of  stochastic  finite  element  methods.  This  is  because  there  is  no 
capability  for  characterizing  or  discretizing  random  fields. 
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FIGURE  2.1b:  External  Concentrated  Loads  on  the  Two-Storey  Braced  Frame 
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FIGURE  2.2c:  Section  Forces  for  a  Tubular  Member 
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FIGURE  2. 2d:  Geometric  Parameters  of  a  Plane  K-joint 
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3.  CASE  STUDIES  OF  SFEM-BASED  COMPUTER  PROGRAMS 

3.1  Introduction 

As  of  the  time  of  writing  this  document,  there  are  only  a  few  computer  programs  available 
in  the  public  domain  for  performing  SFEM-based  reliability  analysis.  The  two  most  prominent  are 
NESSUS  developed  by  Southwest  Research  Institute  (SwRI)  in  San  Antonio,  Texas  on  behalf  of  NASA 
(Lewis  Research  Center)  and  the  CALREL-FEAP  system  developed  at  the  University  of  California  at 
Berkeley.  In  this  chapter,  a  discussion  of  the  main  features,  structure,  and  capabilities  of  both  systems 
are  discussed. 

3.2  The  NESSUS18  System 

3.2.1  General  Description 

The  NESSUS®  software  system  (acronym  for  Numerical  Evaluation  of  Stochastic  Structures 
Under  Stress)  is  a  general  purpose  probabilistic  computer  program  for  predicting  stochastic  structural 
response  due  to  random  loads,  material  properties,  part  geometry,  and  boundary  conditions.  This 
program  is  the  outcome  of  an  ongoing  effort  at  Southwest  Research  Institute,  San  Antonio,  Texas  under 
contract  to  the  National  Aeronautics  and  Space  Administration,  Lewis  Research  Center  (NASA/LeRC). 
The  project,  referred  to  as  "Probabilistic  Structural  Analysis  Methods"  (PS AM)  for  Select  Space 
Propulsion  System  Components  is  being  carried  out  in  two  five-year  phases.  The  first  five-year  PSAM 
program  (1985-1989)  developed  tools  for  predicting  the  probabilistic  structural  response  of  large  scale 
structural  components.  The  second  five-year  program  (1990-1994)  is  extending  this  capability  to 
include  more  general  failure  definitions,  and  addressing  more  complicated  issues  such  as  reliability 
analysis  for  structural  systems,  certification,  and  health  monitoring  [50]. 

The  NESSUS  system  has  several  software  components.  These  include  an  EXPERT  system, 
a  probabilistic  finite  element  code  (PFEM),  a  probabilistic  boundary  element  code  (PBEM),  a  fast 
probability  integrator  (FPI)  module,  and  a  pre-processor  module. 
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The  NESSUS/EXPERT  module  is  an  interactive  menu  driven  expert  system  that  provides 
information  to  assist  in  the  use  of  the  probabilistic  finite  element  code  NESSUS/PFEM  and  the  fast 
probability  integrator  NESSUS/FPI.  This  module  has  the  objective  of  capturing  and  utilizing  PSAM 
knowledge  and  experience  for  guiding  and  assisting  the  engineer  in  the  effective  and  efficient  use  of 
probabilistic  structural  analysis  technology. 

Since  the  PSAM  project  is  designed  to  perform  probabilistic  structural  analysis  on  realistic 
engineering  structures  under  complex  loading  environments,  finite  element  analysis  has  a  key  role  to 
play.  The  NESSUS/FEM  module  is  a  nonlinear  probabilistic  finite  element  code  for  determining 
structural  response  and  sensitivities.  This  program  module  is  based  on  a  mixed  finite  element 
variational  formulation  in  which  the  nodal  degrees  of  freedom  includes  stresses,  strains,  as  well  as 
displacements.  NESSUS/FEM  includes  special  "perturbation"  algorithms  to  obtain  the  response  to  a 
perturbed  set  of  random  variables.  Perturbation  algorithms  provide  efficient  ways  (compared  to  the 
finite  difference  method)  of  computing  the  sensitivity  of  the  structural  response  to  small  fluctuations 
of  the  random  variables  about  a  deterministic  value.  Analysis  types  that  can  be  handled  by  this  program 
module  include  static,  natural  frequency,  buckling,  harmonic  excitation,  and  random  vibration.  The 
element  library  includes  beam,  plate,  plane  stress,  plane  strain,  axisymmetric,  and  3-D  solid  elements. 
Random  variables  may  include  geometry  parameters,  loads,  and  material  properties.  The  program  also 
has  the  capability  to  deal  with  material  and  geometric  nonlinearities. 

NESSUS/FPI  is  a  fast  probability  integration  analysis  program  that  evaluates  point 
probability  estimates  and  cumulative  probability  distributions  (CDF)  for  selected  structural  response 
variables.  The  FPI  procedure  implemented  in  NESSUS  is  based  on  an  extension  of  the  Rackwitz- 
Fiessler  algorithm  [51].  A  procedure  called  the  Advanced  Mean  Value  First  Order  Method  (AMVFO) 
is  used  to  couple  the  finite  element  and  FPI  codes  to  efficiently  obtain  CDF’s  of  large  scale  complex 
engineering  structures  whose  solution  times  are  computationally  expensive.  Improved  versions  of  this 
procedure  (such  as  the  AMV+  method)  have  also  being  developed  and  implemented  into  NESSUS  by 
SwRI  personnel.  The  AMV+  method  is  an  iterative  algorithm  designed  to  improve  the  accuracy  of  the 
advanced  mean  value  method.  This  is  achieved  by  using  the  exact  most  probable  point  locus  (MPPL) 
to  define  the  higher  order  terms  in  the  Taylor  series  expansion  of  the  performance  function  about  the 
mean  value  [50].  Two  different  AMV-t-  algorithms  are  available  in  NESSUS.  The  AMV+  p-level 
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procedure  is  used  to  compute  the  response  corresponding  to  a  specified  probability,  while  the  AMV  + 
z-level  procedure  is  for  the  computation  of  the  probability  corresponding  to  a  prescribed  overall 
response  z(x).  NESSUS/FPI  also  includes  a  Monte-Carlo  simulation  scheme.  Further  details  of  the 
AMV/AMV+  algorithms  are  given  in  Subsection  3.2.2. 

A  schematic  representation  of  the  NESSUS  (Version  5.0)  code  structure  is  shown  in  Figure 
3.1.  In  addition  to  a  library  of  precoded  functions,  user  programmable  subroutines  are  provided  to 
allow  any  performance  function  to  be  defined.  The  user  has  an  option  of  selecting  the  method  of 
structural  response  calculation  (NESSUS/FEM,  user  subroutine,  user  input,  etc.).  One  of  the  new 
features  of  NESSUS  is  the  facility  to  store  probabilistic  results  in  a  database  for  later  post-processing. 

New  features  of  NESSUS  include  systems  reliability  analysis,  adaptive  importance  sampling 
and  probabilistic  fault  tree  analysis.  A  facility  for  probabilistic  fracture  mechanics  analysis  is  also  under 
development  [52]. 

The  NESSUS  system  is  available  from  SwRI  for  a  yearly  license  fee  of  about  $20,000 
(Cdn.).  An  additional  first  time  installation  fee  of  about  $3,600  applies  to  all  licenses.  The  annual 
license  includes  usage,  maintenance,  and  support  of  the  software.  Two  variations  of  the  annual  license 
are: 

i.  Maintenance  and  telephone  support  license,  which  includes  software  usage,  unlimited 
telephone  support,  and  periodic  code  update  and  enhancements;  and 

ii.  Maintenance  and  engineering  support  license,  which  includes  software  usage,  120  hours 
of  SwRI  engineering  (including  telephone)  support  and  periodic  code  updates  and 
enhancements.  This  second  option  is  available  in  the  first  year  only. 

There  was  no  access  to  this  code  during  the  course  of  this  work.  However,  its  theoretical 
foundation  and  some  of  its  applications  are  documented  in  conference  proceedings  and  journal  articles, 
many  of  which  were  thoroughly  studied  by  the  principal  investigator  prior  to  and  during  the  course  of 
the  present  contract. 
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NESSUS  is  a  probabilistic  finite  element  analysis  package  that  permits  an  integrated 
response  analysis  and  reliability  analysis  of  complicated  structures.  The  AMV  algorithm  developed  at 
SwRI  is  unique  to  NESSUS  and  has  the  major  advantage  that  it  efficiently  performs  probabilistic 
analysis  by  a  procedure  that  minimizes  the  number  of  limit  state  function  evaluations  in  comparison 
with  conventional  optimization  schemes.  The  algorithm  has  been  successfully  applied  to  many  practical 
problems  in  the  aerospace  industry.  However,  further  rigorous  applications  are  required  to  fully 
establish  the  robustness  of  the  AMV  algorithm. 

3.2.2  NESSUS  AMV/AMV+  Algorithms 


As  was  pointed  out  in  the  preceding  subsection  the  so-called  AMV/AMV  +  algorithms  were 
developed  at  SwRI  and  are  unique  to  the  NESSUS  system.  In  this  subsection,  these  algorithms  are 
elaborated  upon. 

Consider  a  response  (or  performance)  function  g=g(x),  where  x  is  a  vector  of  N  random 
variables  and  g  is  either  explicit  or  implicit.  The  function  g(x)  is  usually  implicit  when  dealing  with 
complex  structures  where  structural  response  is  evaluated  by  means  of  the  finite  element  method.  The 
AMV  procedure  has  the  fundamental  objective  of  constructing  reasonably  accurate  cumulative 
distribution  functions  (CDFs)  of  g,  namely  Fg,  with  a  very  minimum  number  of  performance  function 
evaluations.  A  description  of  the  iterative  procedure  follows: 

i.  Approximate  g  as  a  linear  function  of  X  via  a  Taylor  series  expansion  about  a  representa¬ 
tive  point,  usually  taken  as  the  mean  value; 

ii.  Based  on  the  approximate  explicit  function  obtained  in  step  (i),  use  reliability  methods  to 
approximate  probabilities  in  selected  points  in  the  sample  space  of  g.  This  first 
approximation  to  Fg  is  called  the  mean-value  first-order  (MVFO)  method,  which  is  usually 
inaccurate; 

iii.  Evaluate  the  function  g  at  each  design  point  to  improve  the  estimate  of  Fg  obtained  in  step 
(ii).  This  is  referred  to  as  the  "first  move"  in  the  AMV  method; 

iv.  a.  Obtain  a  linear  approximation  (Taylor’s  series)  to  g  at  each  of  the  updated  design 

points,  and 

b.  Use  a  fast  probability  integration  (FPI)  scheme  to  obtain  point  probability  estimates 
in  the  "next  move"  to  construct  Ft; 
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v.  Evaluate  the  g  function  at  the  design  points  computed  in  step  (iv);  and 

vi.  Repeat  steps  (iv)  and  (v)  to  improve  the  estimate  of  g  until  desired  accuracy  is  achieved. 

When  the  algorithm  is  terminated  at  step  (iii),  it  is  referred  to  as  the  A  MV  method  (or 
AMVFO,  i.e.  AMV  first-order,  the  first-order  referring  to  the  linear  approximation  in  step  (v)).  For 
a  problem  in  which  M  points  are  selected  for  the  CDF  construction  (i.e.  M  probability  levels)  in  the 
sample  space  of  g,  the  AMVFO  method  requires  a  total  of  J0AMV  function  evaluations  where 

Jo""  =  N  +  M  +  1  .  (31) 

This  number  is  based  on  the  assumption  that  the  gradients  of  the  performance  function  are 
computed  using  the  perturbation  method  as  is  done  in  NESSUS.  When  the  algorithm  includes  steps 
(iv-vi),  the  procedure  has  been  referred  to  as  the  AMV+  algorithm  based  on  the  most  probable  point 
locus  (MPPL).  The  number  of  g-fimction  evaluations  for  the  AMV  +  procedure  in  a  solution  requiring 
L  full  iterations  is  given  by 

=  (N+M+l)  +  LM(N  +2)  .  (3-2) 

A  full  iteration  here  is  defined  as  an  iteration  required  to  update  the  linear  approximation  to  g  at  the 
most  current  design  points,  the  calculation  of  point  probability  estimates  corresponding  to  this 
approximation,  and  "moves"  to  correct  the  g-fimction  values  corresponding  to  the  probability  estimates. 
The  process  is  referred  to  as  the  p-level  version  of  the  AMV -I-  algorithm  because  probability  levels 
are  selected  and  appropriate  response  values  are  determined  for  the  probability  levels. 

The  z-level  AMV+  algorithm,  on  the  other  hand,  is  for  estimating  probabilities 
corresponding  to  given  response  levels,  Zq.  The  steps  involved  may  be  summarized  as  follows: 

i-  Use  the  mean-value  first-order  (MVFO)  method  to  obtain  the  intercept  and  slope  of  the 

cumulative  distribution  function  (CDF)  curve  at  the  50%  probability  level; 

ii.  Let  zt(x)  =  Zq  and  compute  the  MPP  and  CDF  using  MVFO  method; 

iii.  Recompute  Zq  at  the  most  probable  point  (MPP); 
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iv.  Use  corrected  CDF  point  and  information  from  step  (i)  to  fit  a  quadratic  CDF  curve,  and 
then  use  this  curve  to  predict  the  probability  level  P0  for  specified  Z0; 

v.  Set  the  value  of  the  CDF  to  P0  and  use  z,(x)  to  find  the  corresponding  response  value  Zo 
and  the  MPP;  and 

vi.  Use  the  MPP  found  as  a  starting  point  for  repeated  iterations  about  Zq  until  the  probability 
converges;  the  solution  is  expected  to  converge  quickly  because  of  the  quadratic  curve 
fitting  scheme. 

The  main  attractive  feature  of  the  AMV/AMV  +  algorithm  is  the  relatively  small  number 
of  PE  calculations  required  to  establish  the  CDF  of  a  given  performance  function  as  compared  to 
conventional  FORM/SORM  algorithms.  The  algorithm  may,  however,  have  difficulties  in  the  analysis 
of  problems  involving  non-monotonic  or  highly  nonlinear  performance  functions.  Some  thoughts  on 
how  to  exploit  the  nice  features  of  the  AMV  algorithms  in  conjunction  with  schemes  to  improve  the 
robustness  of  the  algorithms  were  generated  during  the  course  of  this  work.  These  ideas  form  the  basis 
of  some  of  the  recommendations  for  future  work  given  in  Chapter  7. 

3.3  The  CALREL-FEAP  System 

3.3.1  General  Description 

The  CALREL-FEAP  system  is  a  combined  package  for  performing  integrated  finite 
element/probabilistic  modelling  of  the  reliability  analysis  of  stochastic  structures.  CALREL  is  a 
general  purpose  structural  reliability  analysis  program  that  has  been  developed  by  Prof.  A.  Der 
Kiureghian  and  his  research  group  at  the  University  of  California,  Berkeley  [53],  FEAP  is  a  finite 
element  analysis  package  developed  by  Prof.  R.L.  Taylor  also  of  the  University  of  California  at 
Berkeley.  The  two  programs  are  connected  via  user  defined  subroutines  for  limit  state  function 
evaluations,  the  gradients  of  the  limit  state  functions,  and  the  possibility  of  user  provided  probability 
distributions.  The  development  of  the  integrated  CALREL-FEAP  system  was  directed  by  Prof.  Der 
Kiureghian  [10]  based  on  the  finite  element  reliability  method  (FERM).  Figure  3.2  illustrates  the 
structure  of  the  CALREL-FEAP  system. 
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3.3.2  CALREL-FEAP  Algorithms 

CALREL  has  four  general  techniques  for  computing  the  failure  probability  corresponding 
to  a  given  performance  function.  These  include  FORM,  SORM,  directional  simulation  (with  exact  or 
approximate  surfaces),  and  Monte-Carlo  simulation.  The  second-order  component  reliability  analysis 
is  performed  by  means  of  either  the  point  fitting  method  or  the  curvature  fitting  method,  or  both,  in 
the  standard  normal  space.  Facilities  are  available  for  the  computation  of  first-order  reliability  bounds 
and  PNET  approximation  for  series  systems,  first-order  reliability  sensitivity  analysis  with  respect  to 
probability  distribution  functions  and  limit-state  function  parameters,  and  Monte-Carlo  simulation  for 
general  systems. 

The  program  has  a  library  of  11  standard  probability  distribution  functions  and  provision 
is  made  for  user-defined  subroutines  for  the  addition  of  other  distribution  functions.  Correlations 
between  the  basic  random  variables  are  permitted. 

The  calculation  of  gradients  in  the  CALREL-FEAP  system  is  based  on  the  adjoint  variable 
method  in  conjunction  with  the  use  of  analytically  derived  partial  derivative  element  matrices  and  load 
vectors.  There  is  provision  for  user-supplied  gradients  in  the  system. 

Current  features  include  the  finite  element  reliability  analysis  of  geometrically  nonlinear 
stochastic  structures  using  the  wide  variety  of  the  finite  elements  that  are  available  in  the  FEAP 
program  which  serves  as  a  subroutine  in  this  system.  The  SORM  procedure  for  this  system  is  based 
on  the  point-fitting  paraboloid  method  developed  by  Der  Kiureghian  et  al.  [54].  Further  details  on 
CALREL-FEAP  algorithms  can  be  found  in  reference  [10]. 
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CALREL  Interface  FEAP 
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FIGURE  3.2:  The  Structure  of  the  CALREL-FEAP  System 
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4-  DESIGN  AND  IMPLEMENTATION  OF  THE  MODEL  SFEM  COMPUTER  PROGRAMS 

4. 1  Basic  Requirements  of  Computer  Programs  for  SFEM-Based  Structural  Reliability  Analysis 

There  are  several  components  of  analysis  capabilities  required  of  any  SFEM-based 
reliability  analysis  software  package.  In  the  present  work,  the  finite  element  reliability  method  (FERM) 
developed  by  Der  Kiureghian  and  his  coworkers  [9,  10],  was  selected  because  of  its  generality  and  the 
other  advantages  alluded  to  in  Chapter  2.  The  specific  components  that  constitute  a  FERM  analysis 
capability  are  highlighted  below. 

One  of  the  first  requirements  is  the  capability  for  discretizing  random  fields  into  random 
variables.  This  is  a  necessary  pre-processing  step  as  the  most  advanced  reliability  analysis  algorithms 
are  formulated  on  the  basis  of  random  variable  models.  A  main  program  is  required  to  drive  the  entire 
probabilistic  finite  element  analysis  process.  This  module  will  accept  input,  control  the  flow  of  the 
various  modules  and  algorithms,  and  control  analysis  output. 

A  finite  element  analysis  (FEA)  platform  is  required  for  the  calculation  of  the  structural 
responses.  In  addition  to  this,  a  module  for  the  computation  of  response  gradients  is  required.  This 
is  a  requirement  of  virtually  all  strategies  for  SFEM-based  reliability  analysis,  except  the  Monte  Carlo 
simulation  scheme. 

The  next  major  module  is  the  reliability  analysis  module,  which  performs  the  calculation 
of  failure  probabilities  and  probabilistic  sensitivity  factors.  This  module  requires  the  capability  of 
transforming  generally  correlated  non-Gaussian  random  variables  to  uncorrelated  standard  normal 
variables.  Then  a  fast  probability  integration  (FPI)  scheme  is  required  to  evaluate  the  integral  that 
defines  the  failure  probability,  Pf,  namely: 

P<  -  fa  -  ■  <41> 

where  fx(x)  is  the  joint  probability  density  function  of  the  random  variables  x.  Also,  the  scope  of  the 
present  work  includes  the  calculation  of  various  probabilistic  sensitivity  measures. 
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Finally,  a  limit  state  function  library  is  required  for  the  definition  of  typical  failure  modes 
of  interest.  However,  in  addition  to  this,  a  provision  could  be  made  for  a  user  to  program  limit  state 

functions  of  interest. 

In  the  following  sections,  the  design  and  computer  implementation  of  the  above 
requirements  are  presented.  The  design  of  the  overall  computer  program  was  done  in  a  way  that 
ensures  easy  enhancements  and  implementation  of  new  strategies  in  the  future.  The  selection  of 
specific  algorithms  was  also  done  with  careful  considerations  of  extension  to  more  advanced  capabilities 

that  are  likely  to  be  of  interest  in  the  future. 

4  2  The  Random  Field  Discretization  Module  (RANFLD) 

Program  RANFLD  is  an  important  pre-processor  for  stochastic  finite  element  analysis 
involving  random  fields.  Two  discretization  procedures  were  implemented  in  RANFLD.  these  are  the 
midpoint  method  and  the  nodal  point  method.  These  techniques  were  selected  because  they  are  the 
most  versatile  as  they  are  applicable  to  all  types  of  random  fields.  It  is  recommended  that  material 
property  fields  (namely:  Young’s  modulus,  Poisson’s  ratio,  and  density)  be  discretized  using  the 
midpoint  method,  while  distributed  random  load  fields  could  employ  the  nodal  point  method.  This 
scheme  is  the  most  compatible  with  the  VAST  FEA  system. 

Stochastic  finite  element  (SFE)  meshes  are  similar  to  conventional  finite  element  meshes 
in  that  they  are  used  to  discretize  a  given  domain.  SFE  meshes  are  recommended  to  be  generally 
coarser  (or  at  least  not  finer)  than  the  mechanical  finite  element  mesh.  RANFLD  is  designed  on  the 
basis  that  the  finite  element  model  generation  program  VASGEN  [55]  will  be  utilized  for  the  generation 
of  stochastic  finite  element  meshes.  Provision  is  made,  however,  for  the  user  to  manually  define  an 
SFE  mesh  in  terms  of  a  conventional  finite  element  mesh  without  a  physical  domain  discretization  of 

the  random  field. 

For  the  midpoint  method,  the  centroids  of  the  stochastic  finite  elements  are  automatically 
computed.  The  user  has  the  options  of  specifying  the  centroidal  point  or  the  preferred  location  for 
computing  the  correlation  matrix.  This  point  may  be  specified  in  terms  of  a  node  number  of  a  FE 
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mesh  or  in  terms  of  the  global  (Cartesian)  coordinates.  Provision  is  made  in  RANFLD  for  the  use  of 
the  same  mesh  for  the  SFEM  and  the  FEM.  Multiple  random  fields  are  permitted  but,  in  the  present 
version,  there  is  the  restriction  of  using  the  same  SFEM  mesh  for  material  property  fields.  The 
material  properties  are  classified  into  four  types  as  follows: 

i.  Type  1:  Young’s  modulus; 

ii.  Type  2:  Poisson’s  ratio; 

iii.  Type  3:  mass  density;  and 

iv.  Type  4:  any  other  type. 

This  classification  is  mainly  for  the  information  of  the  gradient  calculation  programs  and  the  "type" 
of  every  random  variable  generated  by  RANFLD  is  tagged  to  that  variable  throughout  the  analysis. 

The  input  to  RANFLD  consists  of  probabilistic  descriptions  of  all  the  random  fields  and 
the  descriptions  of  the  stochastic  finite  element  meshes  to  be  used  for  discretization.  Each  random  field 
is  defined  via  the  specification  of  the  probability  distribution  function,  mean  function,  variance 
function,  and  the  correlation  function.  Acceptable  probability  distribution  functions  include  any  of  the 
16  distributions  available  in  the  COMPASS  [56]  library  and  also  in  the  FORMREL  module.  Four 
commonly  used  analytical  correlation  function  models  were  implemented  as  a  start-up  library  in 
RANFLD.  It  should  be  noted  that  random  fields  can  also  be  described  in  terms  of  the  variance 
function  (5(t))  and  the  scale  of  fluctuation  ( 8 )  instead  of  the  correlation  function.  However,  the  former 
framework  is  usually  employed  in  conjunction  with  spatial  averaging  techniques.  The  RANFLD 
requirement  is  based  on  practical  information  availability  considerations  and  is  by  no  means  a 
restriction  on  its  capabilities.  This  is  especially  so  as  there  exists  a  well  defined  analytical  relationship 
between  the  variance  function  and  the  correlation  function  [7]. 

The  correlation  function  models  implemented  in  RANFLD  during  the  present  work  include 
the  exponential,  Gaussian,  triangular  and  the  second-order  autoregessive  process  models.  These 
schemes  were  implemented  assuming  the  fields  to  be  homogeneous.  (A  random  field  is  called 
homogeneous  if  all  the  joint  probability  distribution  functions  remain  the  same  when  the  set  of 
locations,  5U  52,  ...,  <5n  is  translated  (but  not  rotated)  in  the  parameter  space.)  A  homogeneous  random 
field  is  the  spatial  counterpart  of  a  stationary  random  process. 
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The  exponential  correlation  function  model  is  associated  with  a  first-order  autoregessive  (or 
Markov)  process  and  has  the  form: 


p(-c)  =  exp(-|t  |/d)  , 


(4.2) 


where  t  =  |i  -  |  is  the  distance  between  points  and  ^ ,  and  d  is  the  correlation  distance.  The 

Gaussian  (or  squared  exponential)  model  has  the  form 

p(x)  =  exp  [-(|x|/d)2]  .  (4'3) 

The  triangular  correlation  function  decreases  linearly  from  1  to  0  as  |  x  |  goes  from  0  to  d, 


(l  -  [x  |/d  ,  |  x  |  5»  d 

pW  =  |p  M>d 


(4.4) 


The  second-order  autoregessive  process  correlation  model  has  the  mathematical  representation. 


PO)  = 


1 


+ 


exp(-|x  |/d) 


(4.5) 


Another  important  capability  that  was  provided  for  the  RANFLD  program  is  the  spectral  decomposition 
of  the  covariance  matrix.  This  procedure  (which  is  used  by  NESSUS)  is  only  applicable  to  Gaussian 
fields.  However,  it  has  the  advantage  in  SFEA  of  being  able  to  reduce  the  number  of  random  variables 
by  truncating  the  probabilistic  correlation  modes.  Two  eigenvalue  orthogonalization  subroutines 
(SPECD1  and  SPECD2)  were  implemented  to  achieve  this. 


Also  in  RANFLD,  names  are  automatically  assigned  to  every  random  variable  generated. 
Apart  from  the  fact  that  the  variable  names  are  used  by  the  reliability  analysis  module,  they  serve  as 
useful  identification  for  the  element  groups  and  property  type  (eg.  E,  v  or  p). 


Figure  4.1  illustrates  the  schematic  flow  of  program  RANFLD.  The  input  file  format  is 
described  in  Appendix  A1  and  is  named  PREFX.RFD,  where,  as  in  other  VAST  analysis  files, 
"PREFX"  is  a  five-character  prefix  prompted  for  at  the  beginning  of  the  program  run.  Necessary  files 
are  opened  and  the  centroids  of  the  stochastic  finite  element  meshes  are  computed  (if  required)  using 
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the  main  subroutines  CENTRD.  This  subroutine,  in  turn,  has  several  subsidiary  subroutines  for 
dealing  with  each  element. 

The  current  capability  is  available  for  seven  elements  in  the  VAST  finite  element  library. 
These  include  the  8-noded  thick/thin  shell  element  (IEC=1),  the  two-noded  general  beam  element 
(IEC=3),  the  three-noded  triangular  plate  element  (IEC=4),  the  four-noded  quadrilateral  shell  element 
(IEC=5),  the  two-noded  bar  element  (IEC  =  8),  the  three-noded  triangular  membrane  element  (IEC=9), 
and  the  eight-noded  isoparametric  membrane  element  (IEC=20). 

The  output  of  program  RANFLD  consists  of  the  following: 

i.  An  echo  of  the  random  field  description; 

ii.  Upper  triangular  portion  of  the  correlation  matrix  printed  row-wise;  and 

iii.  Eigenvalues  and  eigenvectors  of  the  covariance  matrix,  arranged  in  descending  order  of 
eigenvalues. 

Two  input  files  are  produced.  A  formatted  output  file  PREFX.LPD  permits  the  user  to 
view  the  results  of  the  random  field  discretization.  A  binary  file  PREFX.RVB  contains  essentially  the 
same  information  and  is  the  file  utilized  by  the  stochastic  finite  element  program  (STOVAST)  for  the 
definition  of  random  variables  generated  from  random  fields. 

4.3  Algorithms  for  Gradient  Computation 

The  most  important  feature  of  a  stochastic  finite  element  based  reliability  analysis  program 
is  the  calculation  of  the  gradients  of  finite  element  structural  responses  which  appear  in  the  performance 
functions.  This  is  because  algorithms  for  first-order  reliability  methods  (FORM)  and  second-order 
reliability  methods  (SORM)  which  are  efficiently  utilized  for  the  approximate  evaluation  of  failure 
probabilities  rely  on  the  accurate  (repeated)  calculation  of  gradients  of  the  limit  state  functions. 

Several  strategies  are  available  for  the  computation  of  the  gradients  of  finite  element 
responses.  These  include  the  finite  difference  method,  perturbation  method,  and  the  adjoint  method. 
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For  the  present  work,  the  adjoint  method  first  suggested  by  Arora  and  Haug  [57]  and  discussed 
extensively  in  the  monograph  by  Haug  et  al.  [58]  was  the  selected  strategy  for  gradient  computation. 
This  method  is  superior  in  terms  of  both  computational  efficiency  and  accuracy  to  the  perturbation 
method  used  by  Hisada  and  Nakagiri  [14]  and  Liu  et  al.  [15,  31],  the  iterative  perturbation  algorithm 
proposed  by  Dias  and  Nakazawa  [59],  and  of  course  the  finite  difference  method.  The  merit  of  the 
adjoint  method  was  emphasized  in  the  work  of  Liu  and  Der  Kiureghian  [10]  and  was  implemented  m 
the  CALREL-FEAP  system.  The  efficiency  of  the  adjoint  method  was  also  demonstrated  by  Reh  et 
al.  [37],  especially  for  problems  involving  a  large  number  of  random  variables. 


Consider  the  most  general  representation  of  the  limit  state  function,  g,  defined  in  terms  of 
response  variables,  namely:  displacements  (U),  strains  (e)  and  stresses  (or),  and  additionally  some  (finite 
element  or  non-finite  element)  basic  (primitive)  variables  (v).  The  mathematical  representation  of  g 

is  then  given  by: 

,  TT  n  (4.6) 

g  =  g(v,U,e,a) 

or 

g  =  g(v,U(v),e(v,U(v)),  a(v,e(v,U(v)))  (4‘7) 


The  gradient  of  g  with  respect  to  a  basic  random  variable  (v),  Vvg,  is  then  given  by. 


dg 

au 


vu 

V 


Vv€ 

de 


da 


(4.8) 


where  VvU,  Vve,  and  Vv<r  are  the  gradients  (with  respect  to  v)  of  the  displacement,  strain  and  stress 
fields,  respectively.  Using  the  chain  rule  of  differentiation  and  the  respective  functional  dependencies 
of  U,  e  and  a,  it  is  easily  shown  that 


Vvg 


,  it 

dv 


dg  |  +  dg  da  . 

de  dv  U  dcr  3v  6 


dg 

do 


da 

de 


de  j 
dv 
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dg 
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+  dg  de,  +  dg 
de  aU  v  do 


da  |  de  . 
de  lv  oU lv 


Vvu 


(4.9) 


Equation  (4.8)  is  the  most  general  representation  of  the  gradient  of  typical  performance  functions  of 
interest  in  structural  mechanics.  It  may  be  more  compactly  represented  as. 
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where 


and 
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vvg 

+ 

II 

^  vvu , 

3g 

—  1 
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iu  » 


(4.10) 


(4.11) 


(4.12) 


For  the  special  case  of  displacement-type  limit  states  the  functional  representation  of  g  is 
of  the  form: 


g  =  g(v,U(v)) 


and  fv  and  f<j  in  Equations  (4.11)  and  (4.12),  respectively  reduce  to: 


and 


=  ii 

0v  ’ 


(4.13) 


(4.14) 


(4.15) 


Similarly,  stress-type  limit  state  conditions  in  engineering  practice  are  usually  defined  in 
terms  of  a  stress  threshold  (typically  yield  stress  of  the  material)  and  some  other  stress  components. 
The  functional  representation  then  takes  the  form: 


g  =  g(v,  o(v,U(v))  , 


and  fv  and  fu  in  Equations  (4.11)  and  (4.12),  respectively  reduce  to: 

f  =  ii  +  ii.  iSLi 

v  dv  do  dv  u 


(4.16) 


(4.17) 


and 
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'  f  -  dg  do  |  (4.18) 

u  '  3a  dU  v  ' 

It  can  be  appreciated  from  the  above  equations  that  the  calculation  of  the  gradients  of  limit 
state  functions  within  the  framework  of  a  finite  element  analysis  procedure  is  a  very  complicated 
process.  A  major  requirement  of  the  procedure  laid  out  above  is  the  calculation  of  the  gradients  of  the 
displacement  vector  with  respect  to  the  basic  random  variables,  that  is  VvU. 


The  direct  method  of  evaluating  Vvg  is  to  compute  fv,  fy,  VvU  and  then  the  product  fuVvU 
as  required  in  Equations  (4.10),  (4. 1 1)  and  (4. 12).  For  a  linear  static  formulation  of  the  finite  element 
method  (which  is  the  scope  of  this  contract),  the  assembled  equilibrium  equations  appear  in  the  form: 

KU  =  F  ,  (4'19) 

where  K  is  the  assembled  stiffness  matrix  and  F  is  the  assembled  force  vector,  both  random  in  general. 
From  this  equation,  it  follows  that: 

KVU  =  —|  -®U  .  <4-20> 

v  dv  u  dv 

The  direct  computation  of  VvU  from  Equation  (4.20)  is  very  time  consuming.  For  a  structure  with  m 
degrees-of-freedom  and  n  basic  random  variables,  this  approach  requires  n  solutions  of  m  simultaneous 
equations  for  each  limit  state  function  of  interest. 

The  adjoint  method  avoids  the  above  costly  computations  by  solving  the  adjoint  equation: 

XTK  =  fy  ,  (4‘21) 

only  once  for  an  adjoint  displacement  vector  X.  The  product  f0VvU  required  in  Equation  (4. 10)  is  then 
directly  computed  from  the  expression 


XTKVvU  =  ftJVvU  =  XT 


3F  i  3K 


dv 


iu 


dv 


lu 


U 


(4.22) 


in 


which  advantage  of  Equations  (4.20)  and  (4.21)  have  been  taken. 


P152713.PDF  [Page:  71  of  174] 


4.9 


The  evaluation  of  fuVvU  then  requires  the  solution  of  one  equilibrium-type  equation  in 
which  the  stiffness  matrix  involved  is  the  same  stiffness  matrix  for  the  structural  response  in  the  current 
reliability  iteration,  the  computation  of  load  gradients,  and  the  gradients  of  the  stiffness  matrix  with 
respect  to  the  basic  random  variables.  For  this  reason,  two  major  modules  for  gradient  computations 
were  designed  and  developed  in  this  work  for  the  calculation  and  assembly  of  stiffness  and  load 
gradients.  A  summary  of  the  entire  adjoint  method  procedure  for  the  calculation  of  limit  state  function 
gradients  is  schematically  illustrated  in  Figure  4.2. 

The  element  gradient  module  ELEMGD  (similar  to  ELEMS1  in  VAST)  directs  the 
computation  of  the  gradients  of  element  matrices,  with  element-specific  subroutines  implemented  in  a 
submodule  ELGRAD  (similar  to  ELSUB1/ELSUB2  in  VAST).  The  assembly  of  the  element  partial 
derivative  matrices  is  performed  by  a  module  ASSEMG  (similar  to  ASSEM1  in  VAST),  and  the 
solution  of  the  adjoint  problem  is  computed  using  a  module  called  STODIS  (similar  to  DISP1  in 
VAST). 

The  second  gradient  module  is  LOADGD  for  the  calculation  of  load  gradients,  with 
element-specific  subroutines  for  calculating  element  distributed  loads  assembled  in  a  submodule 
LDGRAD. 


The  computation  of  gradients  of  concentrated  nodal  loads  and  the  assembly  of  all  load 
gradients  is  done  within  the  load  gradient  module.  Both  of  these  modules  (i.e.  ELEMGD/ELGRAD 
and  LOADGD/LDGRAD)  are  driven  by  a  subprogram  called  GRADEN.  Subprogram  GRADEN 
includes  a  facility  for  computing  the  required  adjoint  force  and  displacement  vectors  for  displacement 
limit  states,  and  directs  the  calculation  of  gradients  in  the  entire  reliability  analysis  process. 

Three  elements  were  implemented  in  the  gradient  modules.  These  are  the  two-noded 
general  beam  element  (IEC=3),  the  fbur-noded  quadrilateral  shell  element  (IEC=5)  and  the  three- 
noded  triangular  membrane  element  (IEC=9).  The  formulation  of  the  element  and  load  gradients  are 
element-specific.  However,  a  typical  finite  element  has  a  stiffness  matrix  representation  of  the  form: 

K  =  f  B T  D  B  dV  , 

°  J  Qc 


(4.23) 
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where  B  is  the  strain-displacement  matrix,  D  is  the  constitutive  matrix,  V  is  the  volume  and  0,  is  the 
domain  over  which  the  element  is  defined.  The  gradient  with  respect  to  random  material  property  v-, 
can  be  expressed  as: 

El  .  f  BT^BdV,  (4-24> 

dv .  ^ 


since  B  is  generally  independent  of  material  properties. 


For  element  IEC=5,  for  example,  the  (isotropic)  constitutive  matrix,  D,  has  the  form 

1  v  0  0  O' 

v  1  0  0  0 


0 

0 


0 

0 


1  — v 
2 

0 


0 


1  —  v 
2k 


0 

0 


(4.25) 


0  0  0 


0 


1  — v 
2k 


where  E  is  the  Young’s  modulus,  v  is  the  Poisson’s  ratio,  and  k  is  the  shear  displacement  factor.  This 
leads  to  the  following  gradients: 
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The  general  reliability  problem  is  usually  formulated  in  terms  of  a  finite  set  of  basic  random 
variables  X  =  (X„  X2,  ...,  XJ  and  a  limit  state  function  g  =  g(X)  in  which  g  is  the  failure  or 
performance  function.  Failure  is  defined  by  the  event  (g(X)  ^  0},  while  {g(X)  >  0}  identifies  a  safe 
state.  The  probability  of  failure,  Pf,  is  defined  as: 
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Pf  =  Prob{g(X)^0)  = 


/ 

gCCKO 


fx(x)dx 


(4.27) 


where  fx(x)  is  the  multivariate  joint  probability  density  function  of  X. 

The  integral  in  Equation  (4.27)  is  in  general  very  difficult  to  evaluate  and  approximate 
procedures  have  evolved  as  practical  tools  for  efficiently  calculating  the  integral.  Prominent  examples 
in  this  connection  include  the  first-order  reliability  methods  (FORM)  and  the  second-order  reliability 
methods  (SORM)  which  are  well  discussed  in  classical  monographs  on  structural  reliability;  see,  for 
example,  Madsen  et  al.  [60]  and  Melchers  [1].  These  procedures  are  based  on  the  evaluation  of  a 
reliability  index  (/?)  from  which  the  failure  probability  (Pf)  can  then  be  computed.  The  failure 
probability  is  related  to  the  reliability  index  (which  is  the  widely  accepted  probabilistic  measure  of 
safety)  through  the  expression: 

Pfi  =  $(-p)  ,  (4-28) 

where  Pf  is  the  FORM  failure  probability,  and  <$  is  the  standard  normal  cumulative  distribution 
function  (CDF).  The  relationship  in  Equation  (4.28)  is  an  important  one  and  must  be  available  in  any 

reliability  analysis  program.  This  relationship  is  graphically  illustrated  in  Figure  4.4. 

The  scope  of  the  present  work  is  limited  to  FORM  computations.  The  FORM  algorithm 
(FORMREL)  implemented  in  the  reliability  analysis  module  (RELAMM)  follows  the  one  implemented 
in  the  Martec  reliability  analysis  package  COMPASS  [56],  The  main  steps  in  the  FORM-based 
computation  of  /?  are  as  follows: 

i.  Transformation  of  the  vector  of  basic  random  variables  from  the  original  x-space  to  the 
standard  normal  u-space  as  shown  in  Figure  4.5; 

ii.  A  search  (usually  in  the  u-space)  for  the  point  (u*)  on  the  limit  state  surface  g(u)  =  0  that 
has  the  highest  joint  probability  density;  this  point  is  popularly  referred  to  as  the  design 
point,  failure  point,  or  the  most  probable  point  (MPP); 

jji  An  approximation  at  the  MPP  of  the  failure  surface  in  u-space,  and 
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A  computation  of  the  distance  from  the  origin  to  the  MPP  called  0  and  hence  the  failure 
probability. 


The  probability  transformation  (step  (i)  above)  is  presently  accomplished  in  four  substeps. 
First,  there  is  a  transformation  from  the  x-space  to  an  (generally  non-Gaussian)  intermediate  y-space 
of  equivalent  normal  variables  y;  by  equating  their  CDFs  and  PDFs,  that  is: 


*  p-N  -  FGO  , 


y.  -  u, 

<«>  =  f(Xi)  ,  i=l,2,..., 

I  %  J 


(4.29) a 

(4.29) b 


where  $  is  the  standard  normal  PDF,  F(x;)  is  the  marginal  CDF  of  the  basic  random  variable  Xj  and 
f(x;)  is  the  corresponding  PDF.  Equations  (4.29)a  and  (4.29)b  lead  to  the  evaluation  of  the  mean  value 
and  standard  deviation  of  the  resulting  equivalent  normal  distribution;  namely: 


^  =  Xi  -  CTyi  ^  F(Xi)] 


fa^lKx,) 


%  = 


(4.30) a 

(4.30) b 


The  next  substep  then  is  to  transform  the  correlation  matrix  from  the  x-space  to  the  y-space. 
This  correlation  matrix  is  then,  in  turn,  transformed  to  independent  standard  normal  variables  via  an 
orthogonal  transformation: 


r  =  TyTy  , 


(4.31) 


in  which  Ty  is  the  orthogonal  matrix  of  eigenvectors  (R)  of  the  correlation  matrix  in  the  intermediate 
space. 


(Please  note  that  the  superscript  T  in  Equation  (4.31)  denotes  the  matrix  transpose  operator.)  Finally, 
the  normalized  variables  y'  are  transformed  to  the  required  standard  normal  u-space  via  a  linear 
transformation: 
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u  =  [S]"1  [Hf  y'  ,  (4'32) 

where  H  is  the  nxn  matrix  of  scaled  eigenvectors  R  and  [S]  is  the  diagonal  matrix  of  the  eigenvalues. 
The  resultant  transformation  from  the  original  x-space  to  the  u-space  is  usually  denoted  by  the 
transformation  operator  T  such  that: 

U  =  T(X)  .  (4'33) 

The  above  probability  transformation  scheme  (illustrated  schematically  in  Figure  4.5)  has  been  verified 
to  yield  very  accurate  results  in  COMPASS  as  is  also  demonstrated  in  Chapter  5.  A  more  robust 
algorithm  that  has  been  recommended  by  Liu  and  Der  Kiureghian  [27]  will  be  considered  in  future 

work. 

The  search  referred  to  in  step  (ii)  above  is  conducted  by  means  of  the  solution  of  an 
optimization  problem.  The  optimization  problem  pertaining  to  the  calculation  of  the  Hasofer-Lind 
reliability  index  in  the  u-space  may  be  posed  as  follows: 

I  T  n  (4.34)3. 

minimize  D  =  yu*  Uj  =  p  , 

subject  to  g(Uj)  =  0  .  (4.34)b 


The  solution  of  this  problem  locates  the  MPP  and  the  n-dimensional  position  vector  of  this  point,  U“, 
is  given  by: 


U*  =  -cc*p  , 


(4.35) 


where  a  is  the  unit  normal  vector  at  the  MPP,  i.e. 

a*  =  vgan 

|Vg(U  *)  I  ‘ 


(4.36) 


The  optimization  algorithms  implemented  include  the  HL-RF  algorithm  originally  proposed 
by  Hasofer  and  Lind  [61]  and  later  extended  by  Rackwitz  and  Fiessler  [51]  to  include  distribution 
information.  This  is  currently  the  most  widely  used  method  for  solving  the  constrained  optimization 
problem  in  structural  reliability  [62].  The  method  is  based  on  the  recursive  formula. 
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U, 


k+l 


VgT(Uk)Vg(Uk) 


(Vg  T(Uk)Uk  -  g(Uk))  Vg(Uk) 


(4.37) 


Experience  shows  that  for  most  situations  the  HL-RF  algorithm  converges  rapidly. 


Liu  and  Der  Kiureghian  [62]  introduced  a  modified  HL-RF  procedure,  referred  to  hereafter 
as  the  MHL-RF  algorithm,  which  introduces  a  merit  function  m(U)  to  monitor  the  convergence  of  the 
scheme.  The  non-negative  merit  function  proposed  has  the  form: 

m(U)  =  hv  -  Vg(U)|2  +  A  cg(U)2  ,  (4.38) 

2  |Vg(U)|2  2 

where  c  is  a  positive  constant.  In  the  MHL-RF  procedure,  the  new  iteration  point  is  selected  by  a  line 
search  along  the  direction  vector 

4  =  |t7  *  [Vg(Ut)Uk  -  g(Uk)]  Vg(Uk)T  -  Uk  (4.39) 

|Vg(Uk)|2 

until  a  sufficient  decrease  in  m(U)  is  achieved.  The  merit  function  defined  in  Equation  (4.38)  is 
believed  to  be  a  convenient  guide  for  selecting  the  step  size.  Numerical  examples  were  given  in 
reference  [62]  to  illustrate  the  improved  robustness  of  the  MHL-RF  method  over  the  original  HL-RF 
algorithm.  The  MHL-RF  algorithm  was  implemented  in  this  work  because  of  this  robustness  property 
especially  for  finite  element  limit  state  functions.  The  implementation  of  the  algorithm  is  schematically 
illustrated  in  Figure  4.6. 


4.5  Computation  of  Sensitivity  Measures 


Modern  probabilistic  reliability  analysis  programs  must  include  capabilities  for  computing 
sensitivity  measures.  This  is  because  the  identification  of  the  main  sources  of  uncertainty  which  have 
significant  influences  on  the  reliability  of  a  system  is  as  important  as  the  calculation  of  failure 
probabilities.  The  importance  of  probabilistic  sensitivity  estimation  was  demonstrated  in  a  recent  work 
by  the  authors  [63]. 
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Sensitivity  measures  are  utilized  for  two  major  purposes  in  reliability  analysis.  First,  they 
can  be  used  to  identify  the  variables  or  distribution  parameters  (such  as  means,  standard  deviations  and 
correlations)  which  have  major  influences  on  failure  probability.  Second,  they  may  be  used  to  gam 
insight  into  the  complex  behaviour  of  structural  systems.  In  the  first  connection,  especially,  these 
quantitative  measures  provide  guidance  in  the  assessment  of  the  validity  of  reliability  estimates  and  in 
the  definition  of  the  roles  of  the  random  variables  in  subsequent  analysis.  The  first-order  reliability 
method  (FORM)  provides  an  excellent  tool  for  assessing  sensitivities  because  they  can  be  determined 
with  relatively  little  extra  computational  effort,  once  the  failure  probability  has  been  computed. 

Parametric  sensitivity  factors  express  the  sensitivity  of  the  reliability  to  small  changes  in 
statistical  parameters  that  define  the  distribution  of  the  basic  random  variables  X  or  deterministic 
parameters  of  the  performance  function.  The  FORM  expressions  for  parametric  sensitivities  of  the  first- 
order  reliability  index  (/3)  and  failure  probability  (Pf)  are  as  follows  [60,  63]: 

i.  For  a  distribution  parameter  oj„ 

ap  _  J_U*T  da*  =  Au*toT (x*)  (4.40) 

awj  p  6«i  P 

ii.  For  a  performance  function  parameter  p;, 

jLg(u*) 

ap  _  i  utx  au*  _ 
apt  p  ap.  |vg(u*)| 

In  the  above,  u*  is  the  point  in  the  standard  normal  space  (u-space)  on  the  limit  state 
the  highest  joint  probability  density:  the  most  probable  failure  point  (MPP). 


(4.41) 

surface  that  has 
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ill. 


For  any  parameter  q  (oj;  or  p,). 


dq. 


*(-p)  =  , 

0q.  3q. 


(4.42) 


where  4>  is  the  standard  normal  probability  density  function  and  is  the  normal  cumulative  distribution 
function. 


Importance  factors  express  the  relative  importance  of  the  different  sources  of  uncertainty 
associated  with  the  basic  random  variables  that  define  a  problem.  The  scaled  and  normalized  sensitivity 
of  the  reliability  index  with  respect  to  the  design  point  in  the  space  of  the  original  basic  random 
variables  is  mathematically  expressed  as  [9]: 


y(x)  = 


PVPx 

IDVPJ 


Pg(u*)TJu>x 

DaOi*)1^ 


(4.43) 


where  D  is  the  diagonal  matrix  of  standard  deviations  of  X,  V/3  is  the  gradient  of  /J  with  respect  to  x 
at  the  linearization  point,  a(u*)  is  the  unit  vector  normal  to  the  limit  state  surface  at  the  MPP  directed 
towards  the  failure  set,  and  Ju>x  is  the  Jacobian  of  the  probability  transformation  from  the  u-space  to 
the  x-space.  The  importance  factor,  y2  for  each  random  variable  X;  is  defined  as  the  square  of  the 
corresponding  component  of  the  unit  vector  defined  in  the  preceding  equation. 


The  above  sensitivity  measures  were  implemented  in  the  COMPASS  reliability  analysis 
computer  program  developed  at  Martec  Limited  and  were  made  available  to  the  reliability  analysis 
module  of  program  STOVAST:  RELAMM/FORMREL  in  a  submodule  called  SENSCMP. 

4.6  Other  Development  Efforts 

In  a  SFEM-based  structural  reliability  analysis  program,  the  finite  element  analysis  is  driven 
by  the  reliability  algorithm.  For  this  reason,  the  VAST  FEA  program  was  converted  to  a  subroutine 
(VASTM)  so  that  VAST  is  the  FEA  module  of  the  main  SFEM  program  STOVAST.  The  programs 
developed  during  the  course  of  this  work  should,  therefore,  not  be  looked  upon  as  providing  a 
stochastic  finite  element  analysis  capability  for  VAST.  Rather,  a  stochastic  finite  element  analysis 
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program  (STOVAST)  has  been  developed  in  which  VAST  is  used  as  the  finite  element  solver.  This 
point  is  raised  here  to  emphasize  that  although  developments  in  STOVAST  will  be  geared  towards 
compatibility  with  the  VAST  system,  they  will  remain  as  two  separate  programs,  each  with  its  own 
user  manual.  This  situation  will  most  probably  remain  so  until,  hopefully,  stochastic  finite  element 
analysis  becomes  the  accepted  standard  for  FEA. 

In  order  to  make  VAST  serve  as  the  platform  for  FEA  in  this  work,  certain  modifications 
were  required  in  existing  VAST  modules.  First  the  geometry  input  file  had  to  be  modified  to  permit 
the  specification  of  random  element  properties  and  random  loads.  This  triggered  modifications  to  the 
element  subroutines  (ELEMS 1  /ELSUB 1  /ELSUB2)  and  the  load  module  and  its  element  library 

(LOAD1). 

It  was  also  necessary  to  transfer  some  computation  parameters  to  the  various  modules  and 
to  suppress  printout  of  intermediate  (FE)  reliability  iteration  results.  Limit  state  functions  for 
displacements  and  stresses  were  also  made  available  in  program  STOVAST. 

4.7  Input  Requirements  for  Running  Program  STQ.VAST 

Detailed  descriptions  of  the  input  formats  for  running  the  random  field  pre-processing 
program  RANFLD  and  SFEM-based  reliability  analysis  program  (STOVAST)  are  described  in 
Appendix  A1  and  A2,  respectively.  For  STOVAST,  the  input  data  files  required  consists  of  all  VAST 
FEA  input  files  (modified  as  described  in  Appendix  A2)  and  the  main  input  file  PREFX.SFE.  If 
random  fields  have  been  discretized  using  RANFLD,  then  the  random  variables  generated  in  file 
PREFX.RVB  must  also  be  available. 


Loop  Over  all  Random  Fields 
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Input  Definition:  PREFX.RFD 
(Random  Fields  Description) 


l 


Compute  Centroids 


i 


Generate  Correlation  Matrix 


i 


Assign  Names  to  the  Random  Variables 


Perform  Spectral  Decomposition  of  Covariance 
Matrix  (If  Required) 


1 


Generate  Output  Files: 

PREFX.LPD  for  Formatted  Output 
PREFX.RVB  for  Binary  Output  used  by  STOVAST 

End 


FIGURE  4.1:  Illustration  of  the  Flow  of  Program  RANFLD 
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FIGURE  4.2:  Summary  of  The  Adjoint  Method  for  Gradient  Computation 


FIGURE  4.3:  Illustration  of  the  Flow  of  Program  STOVAST 
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FIGURE  4.4:  Relationship  Between  Failure  Probability  and  Reliability  Index 
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FIGURE  4.5:  Illustration  of  the  Transformation  of  Original  Basic  Variables 
(x-space)  to  Standard  Normal  Variables  (u-space) 
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FIGURE  4.6:  Illustration  of  the  MHL-RF  Optimization  Algorithm 
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>N  PROBLEMS  FOR  PROGRAMS  RANFLD  AND  STOVAST 


Consider  a  stochastic  beam  of  length  L  whose  Young’s  modulus  is  best  modelled  as  a 
Gaussian  random  field,  with  a  correlation  function: 


p(Xj,  Xj)  =  exp 


(5.1) 


where  Ax  is  the  distance  between  points  x;  and  xj;  and  X  is  a  dimensionless  correlation  parameter.  The 
random  field  has  a  mean  value  (RMV)  of  2.90E+04  ksi  and  a  standard  deviation  (STD)  of  1.74E+03 
ksi.  It  may  be  recalled  (see  Chapter  4)  that  Equation  (5.1)  corresponds  to  the  exponential  correlation 
function  model  with  a  correlation  distance  d  =  XL.  For  a  value  of  X  =  0.25  used  by  Mahadevan  [36] 
and  a  beam  length  of  192",  this  gives  a  correlation  distance  d  =  CORD  =  48"  for  the  preparation  of 
the  preparation  of  the  PREFX.RFD  file.  Program  VASGEN  [55]  was  used  to  prepare  the  stochastic 
finite  element  mesh  used  for  the  random  field  discretization. 


The  output  of  program  RANFLD  (PREFX.LPD  file)  for  this  problem  is  shown  in  Table 
5.1,  where  it  may  be  noted  that  the  midpoint  method  was  used.  The  variable  names  for  the  four 
random  variables  corresponding  to  the  four  SFEMs  are  given  in  the  output  file  as  YM-G#001,  YM- 
G#002,  YM-G#003,  YM-G#004,  respectively.  These  names  indicate  that  in  the  SFEA  that  may  later 
be  required,  all  deterministic  finite  elements  within  SFEM  no.  1  are  best  grouped  together  as  the  first 
element  group  (hence  the  designations  G#001 ,  etc.).  The  YM  is  a  reminder  that  this  random  variable 
represents  a  Young’s  modulus. 


The  correlation  coefficient  matrix  is  printed  row-wise  as  noted  and  shown  in  Table  5.1. 
These  results  agree,  of  course,  with  those  of  Mahadevan  [36].  Only  the  upper  triangular  portion  is 
printed  because  a  correlation  matrix  is  symmetric.  Further  down  the  table,  eigenvalues  and 
eigenvectors  of  the  covariance  matrix  (not  correlation  matrix)  are  also  printed  in  descending  order  of 
the  magnitude  of  the  eigenvalues. 
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5.2  Random  Field  Discretization  of  a  Stochastic  Plate 


The  second  random  field  discretization  example  is  a  stochastic  rectangular  plate.  The 
Young’s  modulus  of  the  plate  is  modelled  as  a  homogeneous,  lognormal  field  with  a  mean  value  of 
1 .0E+06  units  and  coefficient  of  variation  of  10% .  The  autocorrelation  coefficient  ftmction  is  assumed 


to  have  the  isotropic  form: 


p(Ax)  =  exp 


'Axt  AxN 

,  (0.25L)2 , 


(5.2) 


where  L  is  the  width  of  the  plate  (=16  units).  Equation  (5.2)  represents  the  Gaussian  correlate 
model  with  a  correlation  distance  of  4.0.  This  random  field  is  discretized  using  a  mesh  of  64  stochasttc 
finite  elements  as  shown  in  Figure  5. 1 .  (This  utilized  the  VAST  four-noded  .quadrilateral  shell  element 
(IEC=5).)  The  result  includes  a  correlation  matrix  (64  x  64)  written  into  the  binary  output  file  used 
in  the  stochastic  finite  element  analysis  problem  described  under  Section  5.4. 


5.3  SFF.M  Reliability  Analysis  of  a  Stochastic  Beam 

The  stochastic  beam  shown  in  Figure  5.2  was  first  analyzed  for  the  verification  of  program 
STOVAST.  The  beam  has  a  random  Young’s  modulus,  cross-sectional  area  and  second  moment  of 
area  and  is  subjected  to  a  random  field  load.  The  random  field  has  properties  earlier  described  in 
Section  5.1  above.  The  beam  was  first  analysed  using  a  single  stochastic  finite  element  to  discretize 
the  load  field.  This  is  the  same  as  considering  the  random  field  as  a  random  variable  (WLOAD). 
Table  5.2(a)  gives  a  summary  of  the  basic  variables  involved  in  this  problem.  The  cross-sectional  area 
is  Gaussian  with  a  mean  value  of  7.68  in2  and  a  coefficient  of  variation  (COV)  of  5%,  the  second 
moment  of  area  I.  is  also  normally  distributed  with  a  COV  of  5%.  The  Gaussian  load  has  a  mean 
value  of  0.08  k/in  and  a  COV  of  20%.  The  stochastic  finite  element  reliability  analysis  is  conducted 
with  reference  to  the  limit  state  that  the  vertical  deflection  at  the  free  end  of  the  beam  exceeds  0.5% 
of  the  beam  length.  This  translates  to  a  limit  state  function 
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8  -  ^  -  u/  (5-3) 

where  u^  (=  0.5%  of  L)  is  the  displacement  threshold  and  u/  denotes  the  vertical  displacement  at 
node  5  of  the  structural  finite  element  mesh.  The  results  shown  in  Table  5.2b,  indicate  a  reliability 
index,  /?,  of  1.8627  and  a  failure  probability  of  0.31253E-01.  Importance  factors  are  also  shown  in 
Table  5.2b.  These  answers  are  exact  (as  demonstrated  below)  indicating  the  STOVAST  solution  is 
perfect.  Four  finite  elements  were  used  to  model  the  beam. 


The  closed-form  expression  for  the  tip  deflection  of  a  cantilever  beam  of  length  L  subjected 
to  a  uniform  distributed  load  w  unit  length  is  given  by: 


uy 


WL4 
8  El  ’ 


(5.4) 


where  E  is  the  Young’s  modulus  and  I  is  the  second  moment  of  area.  A  reliability  analysis  was  earlier 
conducted  for  this  beam  using  COMPASS  based  on  the  same  limit  state  condition  defined  in  Equation 
(5.3)  but  with  the  tip  displacement  computed  using  the  explicit  representation  in  Equation  (5.4).  The 
results  are  exactly  the  same  as  those  summarized  above  and  detailed  in  Table  5.2(b)  for  both  /S  and  the 
probabilistic  sensitivity  factors.  A  detailed  STOVAST  output  for  the  4-SFEM  case  is  shown  in 
Appendix  Bl. 


This  STOVAST  result  not  only  illustrated  the  successful  solution  of  a  stochastic  finite 
element  problem  but  demonstrates  two  important  points  in  the  development  and  application  of  SFEM 
procedures.  The  first  point  is  that  if  accurate  (deterministic)  FEA  solutions  can  be  found  (as  was  the 
case  for  this  beam  problem),  then  for  well-behaved  limit  state  functions  such  as  the  present  one,  very 
accurate  reliability  estimates  can  be  expected.  The  other  notable  point  is  the  power  and  utility  of 
response  surface  methods.  Equation  (5.4)  is  a  response  surface  representation  of  the  finite  element 
random  variable  u5y  -  i.e.  the  vertical  displacement  of  node  5.  Thus,  for  more  complicated  models, 
if  reliable  (i.e.  reasonably  accurate)  response  surfaces  can  be  constructed  in  the  neighbourhood  of 
critical  structural  locations,  then  reliability  analysis  results  obtained  can  be  depended  upon. 
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5.4  SFEM  Reliability  Analysis  of  a  Stochastic  Plate 

The  stochastic  plate  described  in  Section  5.2  was  also  analyzed  for  STOVAST  verification. 
The  plate  with  a  lognormal  homogeneous  random  field  Young’s  modulus  is  subjected  to  two  random 
point  loads  in  the  global  x  and  y  directions  [64],  The  loads  have  mean  values:  and  p,  of  4,000  and 
1,000  units,  respectively  and  COVs  of  10%  each.  Both  loads  are  assumed  to  be  normally  distributed 
but  statistically  independent.  The  Poisson  ratio  is  assumed  to  be  deterministic  (i.e.  fixed)  with  a  value 
of  0.2.  The  failure  criterion  examined  was  the  exceedance  of  the  horizontal  displacement  at  the  point 
of  load  application  above  a  threshold  of  0.03  units. 

Three  analyses  were  conducted  for  this  problem,  all  utilizing  a  mechanical  FE  mesh  of  64 
quadrilateral  shell  elements.  The  first  analysis  considered  the  Young’s  modulus  random  field  as  a 
single  random  variable  and  yielded  a  reliability  index,  ;5,  of  1.7714  with  a  corresponding  failure 
probability  of  3.825%.  A  detailed  printout  of  STOVAST  output  file,  PREFX.LPT  for  this  analysis 
is  provided  in  Appendix  B2.  This  output  contains  the  description  of  the  input  random  variables,  finite 
element  input  data,  intermediate  FE  reliability  iteration  results  for  both  structural  response,  loads,  and 
response  gradients,  and  the  reliability  analysis  results.  (Please  note  that  some  of  the  output  data  were 
printed  out  to  monitor  the  convergence  performance  of  the  FPI  algorithms  and  a  provision  would  be 
made  in  future  (i.e.  release)  versions  to  suppress  some  of  the  output  information.) 

The  second  analysis  was  conducted  to  check  the  robustness  of  the  probability  transformation 
scheme  for  highly  correlated  random  variables.  To  achieve  this,  the  random  field  was  discretized  using 
64  SFEMs  but  assuming  perfect  correlation  structure  (i.e.  pX}  =  1.0)  between  the  resulting  random 
variables.  This  correlation  option  is  the  null  correlation  choice  in  the  RANFLD  correlation  function 
library  and  may  be  selected  by  setting  the  parameter  ICOREM  (see  Appendix  Al)  to  zero.  This 
analysis  gave  exactly  the  same  reliability  results  (as  it  should)  as  obtained  in  the  first  analysis, 
demonstrating  that  RELAMM  can  adequately  handle  very  closely  correlated  (in  this  case  perfectly 

correlated)  random  variables. 

The  third  analysis  conducted  for  this  problem  used  a  64  SFEM  mesh  to  discretize  the 
Young’s  modulus  random  field  but  this  time  utilizing  the  actual  correlation  function  model  described 
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in  Section  5.2,  that  is  Equation  (5.2).  This  analysis  yielded  a  reliability  index  of  /3  =  2.1139  or  Pf  = 
1.726%.  We  note,  as  expected,  that  the  value  of  the  reliability  index  has  gone  up  with  reference  to 
the  other  analyses.  This  is  because  the  single  random  variable  model  of  the  random  field  over¬ 
represents  the  variability  in  the  parameter  and  hence  a  higher  failure  probability.  This  last  analysis  also 
demonstrates  that  random  variable  models  of  parameters  that  are  truly  random  fields  give  very 
conservative  estimates  of  safety.  The  example  also  underlines  the  importance  of  random  field 
discretization  capabilities  and  the  merit  of  the  SFEM  approach  itself.  All  the  above  problems  were 
solved  using  both  the  HL-RF  and  MHL-RF  algorithms. 

5.5  SFEM  Reliability  Analysis  of  a  Perforated  Stochastic  Plate 

Figure  5.3  shows  a  quarter  finite  element  model  of  a  stochastic  plate  with  a  circular  hole. 
The  square  plate  is  of  dimensions  24"  x  24"  with  a  central  hole  of  radius  2".  The  plate  is  subjected 
to  a  random  distributed  load  at  the  edges.  This  load  is  assumed  to  be  normally  distributed  with  a  mean 
value  of  (2  ksi  and  a  cov  of  20%).  The  Young’s  modulus  of  the  plate  is  considered  to  be  a 
homogeneous  random  field  with  a  mean  value  of  0.290E+05  ksi  and  a  cov  of  10% .  The  correlation 
of  the  field  is  assumed  to  be  Gaussian,  with  a  correlation  distance  of  6". 

The  discretization  (mechanical  and  stochastic)  was  performed  with  the  four-noded 
quadrilateral  shell  element  as  shown  in  Figure  5.3.  The  four  bold  elements  denote  the  stochastic  finite 
elements  used  for  random  field  discretization,  while  the  24  lighter  elements  denote  the  mechanical  finite 
elements  used  for  discretizing  the  continuum  structure. 


Stochastic  finite  element  reliability  analysis  of  the  plate  is  conducted  with  reference  to  the 
displacement  limit  state  function: 

g(x)  =  nthr  ~  us  (5-5> 

where  uUir  is  the  displacement  threshold  of  0.600E-02"  and  u/  denotes  the  horizontal  displacement  at 
finite  element  node  number  20  (see  Figure  5.3)  of  the  plate.  The  reliability  analysis  results  are: 
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(3  =  2.5746 

pf  =  0.501 83 E-02 

The  importance  factors  are  27.71  %  for  the  combined  random  variables  describing  the  Young’s  modulus 
and  72.29%  for  the  random  load. 
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TABLE  5.1:  RANFLD  Output  for  Stochastic  Beam  Problem 

*  *  *  RANDOM  FIELD  DISCRETIZATION  *  * 

PREFIX  DEFINING  SFEM-REL.  ANALYSIS  FILES  r  WALS8 
STOCHASTIC  BEAM  WITH  RANDOM  YOUNGS  MODULUS 
TOTAL  NO.  OF  RANDOM  FIELDS  -  1 

NO.  OF  SFEM  MESHES  USED 

***  DISCRETIZATION  OF  RANDOM  FIELD  NO.  1 

ICATFL  -  1  ITYPE  -  1 

RANDOM  VARIABLE  CLASSIFICATION  :  KRVT  -  2 

ICDF  -  2 

RMV  -  0.290000E+05  STD  -  0.174000E+04 

ICOREM  -  1  IRHOSP  -  I  CORD  -  0.48CS-02 

PREFIX  DEFINING  SFEM  MESH  FOR  THIS  FIELD  :  WALS8 
CENTROIDAL  COORDS.  OF  SFEMs  COMPUTED  AUTOMATICALLY 
CORRELATION  COEFFICIENT  MATRIX  PRINTED  ROWWISE 
NOTE  :  ONLY  UPPER  TRIANGULAR  PORTION  IS  PRINTED 

VARIABLE  NUMBER  :  1 

VARIABLE  NAME  :  YM-G#OOI 
ROW  #  1  HAS  4 

0 . 100000E+01  0.367879E+00 

VARIABLE  NUMBER  :  2 

VARIABLE  NAME  :  YM-G#002 
ROW  #  2  HAS  3 

0 . 10Q000E+01  0.367879E+00 

VARIABLE  NUMBER  :  3 

VARIABLE  NAME  :  YM-G#003 
ROW  #  3  HAS  2 

0 . lOOOOOE+Ol  0 .367879E+00 
VARIABLE  NUMBER  :  4 

VARIABLE  NAME  :  YM-G#004 
ROW  #  4  HAS  1 

0 . 100000E+01 

EIGENVALUES  AND  EIGENVECTORS  OF  COVARIANCE  MATRIX  : 

***  PROBABILISTIC  MODE  NO.  «**  :  I 

EIGENVALUE  -  0.525768E+07 

CORRESPONDING  EIGENVECTOR  : 

0.417922E+00 
0.570387E+00 
0 . 570387E+00 
0 . 417922E+00 

***  PROBABILISTIC  MODE  NO.  ***  :  2 

EIGENVALUE  -  0 . 324831E+07 

CORRESPONDING  EIGENVECTOR  : 

-0.625406E+00 
-0.329950E+00 
0. 329950E+00 
0.62S406E+00 

***  PROBABILISTIC  MODE  NO .  *  *  *  :  3 

EIGENVALUE  -  0.206205E+07 

CORRESPONDING  EIGENVECTOR  : 

-0.570387E+00 
0.4 17922E+00 
0 . 4 17922E+00 
-0 . 570387E+00 

***  PROBABILISTIC  MODE  NO.  ***  ;  4 

EIGENVALUE  -  0.154237E+07 

CORRESPONDING  EIGENVECTOR  : 

0 . 329950E+00 
-0. 62S406E+00 
0. 625406E+00 
-0 .32995OE+O0 


MEMBERS 

0 . 135335E+00  0.497871E-01 


MEMBERS 
0 . 135335E+00 


MEMBERS 


MEMBERS 
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TABLE  5.2a:  Random  Variables  Describing  the  Stochastic  Beam  Problem 
(RELAMM  is  the  Reliability  Analysis  Module  of  STOVAST) 


— 

— 

— 

RELAMM 

VARIABLE  LISTING 

NO. 

|  NAME 

i 

DISTR . TYPE  1 

MEAN  1 

STAND .DEV. 

|  DISTR.PAR. 

i 

|  TMOD 

i 

NORMAL  1 

0 . 2900E+05  | 

| 

0.1740E+04 

|  0.2900E+05 
|  0.1740E+04 

2 

1 

|  AREA 

1 

NORMAL  i 

0.7680E+01  | 

j 

0  .3840E+00 

|  0.7680E+01 
i  0.3840E+00 

3 

1 

|  AIZ 

1 

t 

NORMAL  1 

0.3010E+03  | 

| 

0 . 1505E+02 

|  0 .3010E+03 
|  0.1505E+02 

4 

5 

6 

7 

i 

|  PR 

|  DEN 

|  ZERO 

|  WLOAD 

i 

i 

i 

i 

i 

FIXED  1 
FIXED  I 
FIXED  1 
NORMAL 

0.3000E+00  | 
0.1000E+01  1 
O.OOOOE+OO  1 
0.8000E-01  | 

| 

0.1600E-01 

|  0.3000E+00 
[  0.1000E+01 
|  O.OOOOE+OO 
|  0 . 8000E-01 
|  0.1600E-01 

8 

1 

|  UTHR 

1 

1 

FIXED  I 

0.9600E+00  | 

|  0.9600E+00 
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STOCHASTIC  FINITE  ELEMENT  DISCRETIZATION  TESTING 


81  STRUCTURAL 


FINITE  ELEMENT 

MODEL 

ELEMENT  GROUPS 

ALL 

2 

X 

1 _ 

4 . 978E+00 


FIGURE  5.1:  Random  Field  Discretization  of  a  Stochastic  Rectangular  Plate 
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FIGURE  5.2:  Structural  and  Stochastic  Finite  Element  Meshes  for  a  Stochastic  Beam 
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MECHANICAL  FINITE  ELEMENT 


STOCHASTIC  FINITE  ELEMENT 


FIGURE  5.3:  Mechanical  and  Stochastic  Finite  Element  Meshes  of  a  Stochastic  Plate 
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6-  BIAS  AND  MODELLING  ERRORS  IN  SFEM-BASED  RELIABILITY  ANALYSTS 

6.1  Classification  of  Uncertainties 

There  are  several  possible  sources  of  uncertainty  when  dealing  with  stochastic  structures 
for  which  we  have  highly  recommended  the  application  of  stochastic  finite-element  methodologies. 
These  include  inherent  randomness,  estimation  errors,  model  imperfections  and  human  error.  Inherent 
randomness  is  usually  due  to  natural  variabilities  in  physical  phenomena  or  processes.  Examples  of 
these  include  environmental  loads  imparted  on  marine  structures  by  the  inherently  random  ocean 
environment,  material  properties,  member  sizes  and  geometry.  This  class  of  uncertainties  is  usually 
referred  to  as  Category  I  or  aleatory  uncertainties.  The  other  class  of  uncertainties  called  Category  II 
or  epistemic  uncertainties  encompasses  estimation  errors,  model  imperfections,  and  human  errors. 
They  are  generally  due  to  lack  of  sufficient  data,  knowledge  or  understanding  of  a  complex 
phenomenon  (i.e.  ignorance),  or  use  of  simplified  models.  Examples  of  this  latter  category  include 
errors  in  analysis  results  due  to  inadequate  analytical  procedures  for  load/response  computation  and 
human  errors  induced  during  finite  element  model  generation  or  other  input  descriptions. 

There  exists  a  fundamental  distinction  between  aleatory  and  epistemic  uncertainties. 
Whereas  inherent  randomness  is  intrinsic  to  nature  and  usually  beyond  the  engineer’s  control,  epistemic 
uncertainties  are  extrinsic  and  to  some  extent  reducible.  Reducibility  in  the  latter  class  stems  from  the 
fact  that  an  analyst  may  choose  to  obtain  additional  information  to  improve  the  accuracy  of  the 
estimation  procedures.  This  may  be  in  the  form  of  gathering  more  data  for  more  accurate  probabilistic 
characterization  of  governing  engineering  parameters,  the  execution  of  appropriate  experiments  to  gain 
a  better  understanding  of  the  relevant  phenomena,  or  the  development  and  application  of  more 
sophisticated  models  for  analysis.  Uncertainty  due  to  human  error  may  also  be  reduced  by 
implementing  rigorous  quality  control  measures  in  the  analysis  process  [65].  It  is  pertinent  to  point 
out  that  the  avenues  for  reducing  epistemic  uncertainties  are  associated  with  costs  in  terms  of  money 
and  time.  However,  epistemic  uncertainties  are  information  sensitive,  unlike  aleatory  uncertainties. 
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6.2  Quantification  of  Uncertainties 

Two  important  measures  of  uncertainties  are  the  bias  factor  and  the  coefficient  of  variation 
(COV).  The  bias  factor  for  a  parameter  x  with  a  true  value  x*  is  defined  as: 

B  =  ft  (6.1) 

xp  ’ 

where  xp  is  the  predicted  value  of  the  parameter.  This  parameter  measures  the  systematic  error  in  the 
prediction  of  the  average  value  of  the  desired  quantity.  The  bias  factor  is  often  applied  constantly  to 
predicted  values  and  may  be  considered  to  be  a  random  variable.  The  COV,  5,  is  the  ratio  of  the 
standard  deviation,  <x„  to  the  mean  value,  I,  and  relates  to  the  standard  error  pertaining  to  the 

statistical  imperfection  in  the  estimated  control  value  (i.e.  mean  or  median): 

x  .  (6-2) 


The  credibility  of  uncertainty  measures  is  important  for  the  obvious  reason  that  the  validity 
of  calculated  structural  responses  or  failure  probabilities  depend  on  these  measures. 


6.3  Analysis  of  Uncertainties 


The  bias  factor  corrects  predicted  results  to  give: 


=  Bx_ 


(6-3) 


xt  — p 

This  model  is  popularly  referred  to  as  the  Ang-Cornell  model  [66],  Ditlevsen’s  model  [67]  has  a  more 
general  representation  of  the  form: 

xt  =  BXp  +  b  <6-4> 


where  b  is  a  Gaussian  random  variable  which  is  statistically  independent  of  x.  The  Ang-Cornell  model 
is,  nevertheless,  often  preferred  in  the  propagation  of  uncertainties  because  it  is  simpler  and  requires 
less  information  in  order  to  determine  the  statistics  of  its  parameters.  The  Ang-Cornell  model  is 
particularly  suitable  for  lognormal  random  variables,  especially  when  expressed  as  a  product  of  several 
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other  variables.  Detailed  illustrations  and  applications  of  this  application  can  be  found  in  the  works 
of  Nikolaidis  and  Kaplan  [68]  and  Bea  [69].  However,  the  basic  principle  for  propagating  uncertainties 
using  the  first-order  second  moment  method  is  summarized  below  following  the  work  of  Bea  [69], 

For  the  addition  or  subtraction  of  two  normal  random  variables,  (X  ±  Y)  =  Z,  the  mean 
of  the  resultant  distribution  can  be  calculated  as: 

Z  =  X  ±  Y  (6.5) 

The  standard  deviation  of  the  resultant  distribution  can  be  calculated  as: 

az  =  \la\+o\±2poxaY  (6-6) 

For  the  multiplication  of  two  random  variables,  (X  Y)  =  Z,  the  mean  of  the  resultant 
distribution  can  be  calculated  as: 


Z  =  X  Y  (6-7) 

The  standard  deviation  of  the  resultant  distribution  can  be  calculated  as: 


az  -  X  Y  V(1  +  P2)  [v**V*  *  (V*V*)]  (6'8) 

for  p  =  0 

vz  -  /v*  +  v-J  +  (vX)  (6'9) 

For  the  division  of  two  random  variables,  (  X  /  Y  )  =  Z,  the  mean  of  the  resultant 
distribution  can  be  calculated  as: 


Z  = 


X 

Y 


(6.10) 


The  standard  deviation  of  the  resultant  distribution  can  be  calculated  as: 
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az  =  r/Y2)  M  +  V*  "  2P (VXVY) 


(6.11) 


or  for  p  =  0 


Vz  = 


Yx  +  4 


(6.12) 


are 


When  the  random  variables  X  and  Y  can  be  considered  independent  (p-0),  and  Vx  and  VY 
small  (V  <  <  1),  then  for  either  the  multiplication  or  quotient  of  two  random  variables: 

(6.13) 


vz  -  /vT+  v* 


For  the  square  of  X  (Z  —  X2): 


Z  =  +  a 


(6.14) 


az  =  X2  az  +  2  oz 


More  generally,  if  Z  =  X  Y“ 


Z  =  X  Y“ 


(6.15) 


oz  =  \ja2x  +  aaY' 


For  lognormal  distributed  variables  (logarithms  of  the  variables  are  normally  distributed): 


ctx  =  and  X  =  X50  exp  (0.5  oj 


(6.16) 


and  relations  similar  to  the  above  can  be  utilized  for  propagating  uncertainties. 

A  framework  for  incorporating  model  uncertainties  in  advanced  reliability  methods 
employing  full  distributional  information  was  developed  by  Der  Kiureghian  [70].  This  work, 
formulated  on  the  Bayesian  approach,  derived  simple  formulas  that  show  the  effect  of  model 
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uncertainty  on  the  reliability  index.  The  Der  Kiureghian  model  is  more  general  than  both  the  Ang- 
Cornell  and  Ditlevsen  models  and  permits  the  calculation  of  the  variance  of  the  reliability  index  to 
reflect  the  influence  of  modelling  errors  or  uncertainties. 

In  our  opinion,  the  Der  Kiureghian  model  is  the  most  realistic  vehicle  for  accounting  for 
bias  and  modelling  uncertainties  in  SFEM-based  structural  reliability  analysis.  This  is  for  the  obvious 
reason  that  stochastic  finite  element  analysis  is  a  sophisticated  process  that  is  difficult  to  express  in 
simple  analytical  forms  such  as  are  required  by  the  Ang-Cornell  and  Ditlevsen  models. 

6-4  Sources  of  Bias  and  Modelling  Errors 


The  sources  of  bias  and  modelling  errors  in  SFEM-based  structural  reliability  analysis 
include  the  following: 


Data  collection/characterization  for  material  properties  and  structural  geometry; 

Probabilistic  description  of  random  fields,  i.e.  whether  the  most  suitable  probability 
distribution  function  or  correlation  function  has  been  used  to  define  the  field.  Errors  in 
these  can  be  reduced  by  collecting  more  data  and  by  validation  of  the  stochastic  character 
of  the  field.  Validation  can  be  accomplished  by  the  reverse  process  of  generating  Monte 
Carlo  simulations  of  the  random  field  and  comparing  these  to  measured  values  of  this  field 
at  selected  locations; 

Discretization  of  random  fields.  Some  probabilistic  information  is  lost  when  the  field  over 
a  domain  is  represented  by  pointwise  values  (i.e.  random  variables); 

Probability  transformation  scheme.  The  transformation  of  generally  correlated  variables 
to  standard  normal  (statistically  independent  variables)  involve  approximations,  and  hence 
some  information  is  also  lost  in  this  process; 

Finite  element  model  generation.  Discretization  errors  are  incurred  in  the  process  of 
converting  a  continuum  structure  into  several  discrete  elements; 

Finite  element  formulation.  These  are  errors  due  to  inappropriate  use  of  finite  elements 
or  boundary  conditions.  This  may  also  be  due  to  inadequate  modelling  of  the  physical 
phenomena,  for  example,  using  linear  instead  of  nonlinear  finite  element  analysis,  failing 
to  utilize  inappropriate  constitutive  models,  or  failing  to  account  for  interaction  of  certain 
phenomena  (eg.  thermomechanical  coupling  in  processes  involving  large  inelastic 
deformations); 
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vii.  Calculation  of  FE  structural  response  such  as  displacements  and  stresses.  These  are  the 
solution  errors  in  FEA  which  may  sometimes  be  significant  (eg.  in  nonlinear  FE  analysis); 

viii.  Calculation  of  gradients  of  the  structural  response; 

jx  Calculation  of  the  failure  probability  (and  hence  the  reliability).  These  calculations  usually 

employ  approximate  procedures  which  introduce  errors  in  the  computed  reliability  indices, 

x  Mathematical  definition  of  the  failure  conditions  (i.e.  limit  states)  of  interest, 

xi.  Definition  of  loads;  and 

xii.  Human  error. 

For  marine  structures,  Nikolaidis  and  Kaplan  conducted  [68]  a  study  dealing  with 
uncertainties  in  stress  analysis  on  behalf  of  the  U.S.  Ship  Structure  Committee.  The  study  dealt  mainly 
with  load  characterization  and  response  computation,  but  did  not  deal  with  reliability  estimation.  The 
above  list,  although  generic  in  form,  is  therefore,  more  comprehensive  as  the  scope  of  stochastic  finite 
element  analysis  is  much  broader  than  the  focus  of  the  study  in  reference  [68]. 

As  we  pointed  out  in  the  preceding  subsection,  the  work  of  Der  Kiureghian  [70]  is 
recommended  as  a  framework  for  realistically  accounting  for  model  uncertainties  in  SFEA.  Another 
strategy  that  is  worthy  of  consideration  is  the  use  of  confidence  bounds  on  the  computed  values  of 
failure  probabilities  (or  reliability  indices).  This  is  clearly  an  area  that  would  benefit  from  further 
research  investigations. 
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7-  CONCLUSIONS  AND  RECOMMENDATIONS 

A  research  and  development  program  for  the  application  of  stochastic  finite  element-based 
reliability  analysis  has  been  successfully  executed  in  this  work.  The  detailed  state-of-the-art  literature 
review  of  various  SFEM  methodologies  permitted  a  solid  background  for  the  evaluation  of  the  different 
strategies  and  substrategies  pertaining  to  the  diverse  technical  requirements  of  a  probabilistic  finite 
element  analysis.  Case  studies  of  the  CALREL-FEAP  system  and  the  NESSUS  system  permitted 
practical  evaluation  of  some  of  the  available  strategies. 

The  random  field  modelling  capabilities  implemented  in  RANFLD  are  the  most  versatile 
because  they  are  applicable  to  non-Gaussian  fields  and  have  no  restrictions  on  the  geometrical  shape 
of  the  stochastic  finite  element  meshes.  A  respectable  library  of  correlation  models  were  made 
available  to  permit  easy  applications.  Furthermore,  the  RANFLD  strategy  does  not  suffer  from 
numerical  stability  problems  that  may  result  in  the  reliability  analysis  module  because  the  correlation 
coefficient  matrix  is  always  positive  definite. 

The  finite  element  reliability  method  (FERM)  implemented  in  program  STOVAST  under 
this  contract  was  demonstrated  to  efficiently  give  very  accurate  results.  The  verification  problems 
solved  in  Chapter  5  demonstrate  not  only  the  overall  proper  operation  of  program  STOVAST  but  its 
various  algorithms  such  as  those  for  gradient  computations,  fast  probability  integration,  and  sensitivity 
computations.  Also  very  important,  program  STOVAST  was  designed  to  permit  easy  extension  of  its 
capabilities  either  in  terms  of  additional  reliability  analysis  procedures  or  the  implementation  of  new 
SFEM  strategies. 

The  sources  of  bias/modelling  errors  and  uncertainties  in  the  application  of  SFEM-based 
structural  reliability  analysis  were  also  examined  during  this  work.  Available  procedures  for 
propagating  uncertainties  in  engineering  analysis  were  briefly  reviewed  and  some  suggestions  were 
made  on  a  suitable  approach  for  the  SFEM  process. 

The  verification  of  program  STOVAST  and  the  research  studies  conducted  in  this  work 
have  produced  very  encouraging  results.  Based  on  this  observation,  and  noting  the  unique  promise  of 
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the  stochastic  finite  element  method  as  the  most  viable  strategy  for  realistically  accounting  for  the 
various  uncertainties  present  in  ship  structural  analysis,  it  is  recommended  that  further  work  be 

pursued. 

First,  it  is  desirable  to  consolidate,  expand,  and  further  verify  the  existing  capabilities  of 
program  STOVAST.  This  will  include  the  addition  of  more  elements  to  the  stochastic  finite  element 
library  and  the  upgrading  of  the  routines  in  both  RANFLD  and  STOVAST  for  compatibility  with  the 
nonlinear  version  of  VAST  (Version  7).  The  capability  for  modelling  random  nodal  coordinates  was 
not  implemented  because  of  the  additional  complications  in  the  calculation  of  the  partial  derivative  of 
load  vectors  and  element  matrices.  This  capability  should  be  provided  in  future  work  so  that  a  full 
stochastic  structural  geometry  representation  of  the  mechanical  finite  elements  will  be  permitted.  A 
more  fully  automated  facility  for  the  creation  of  modified  (STOVAST)  PREFX.GOM  files  from  results 
of  random  field  discretization  and  the  conventional  VAST  geometry  files  will  also  be  useful. 

Further  enhancements  to  the  reliability  analysis  algorithms  would  also  be  beneficial.  We 
believe  the  development  of  the  Nataf  model  for  the  probability  transformation  from  x-space  to  u-space 
will  provide  a  more  robust  platform.  An  additional  optimization  algorithm  suitable  for  locating  the 
MPP  for  large  scale  problems  should  be  implemented.  This  will  further  complement  the  existing  HL- 
RF  and  MHL-RF  algorithms.  Also,  information  on  probabilistic  sensitivity  factors  could  be  used  to 

speed  up  the  reliability  iterations. 

The  second  area  that  would  require  further  work  is  the  development  of  a  new  (alternative) 
SFEM  strategy.  Some  advanced  algorithmic  concepts  which  have  been  developed,  tested,  and  (in  some 
cases)  commercialized  by  other  leading  researchers  should  be  seriously  considered  for  implementation 
in  STOVAST. 

Recent  and  current  DREA  research  efforts  in  the  area  of  probabilistic  mechanics  in  general 
and  in  the  subject  of  SFEM-based  reliability  analysis  in  particular  represent  novel  attempts  to  apply  this 
new  technology  to  ships.  However,  several  technical  issues  that  are  unique  to  the  size  scale, 
operational  environment,  and  specialized  architecture  of  ship  structures  need  to  be  researched. 
Included  in  this  last  are  topics  such  as  various  procedures  for  the  optimization  of  computational 
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efficiency,  efficient  implementation  of  ship-specific  engineering  failure  modes  (both  serviceability  and 
ultimate)  and  development  of  time-dependent  reliability  models. 

Overall,  we  believe  the  execution  of  this  project  is  a  major  step  towards  the  practical 
application  of  SFEM-based  reliability  analysis  to  ship  structural  integrity  assessment.  It  is  hoped  that 
the  scientific  authority  continues  to  maintain  interest  and  give  the  dedicated  commitment  and  support 
for  the  vigorous  pursuit  of  this  laudable  objective. 
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INPUT  FORMAT  FOR  THE  RANDOM  FIELD  DISCRETIZATION  PROGRAM:  RANFLD 
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Al-1 

FORMAT  OF  RANFLD  INPUT  DATA:  PREFX.RFD 

Card  1  (A72) 

TITLE  =  Title  describing  the  random  field  problem 

Card  2  (15) 

NRFLDS  =  Number  of  random  fields  involved  in  the  problem 

Card  3  (15) 

NMESH  =  Number  of  stochastic  finite  element  meshes  used 

|  Cards  4-11  are  to  be  provided  NRFLDS  times. 

Card  4  (A8) 

RANFIELD  =  Eight  character  literal  constant  as  header  for  the  beginning  of  a  random 
field  description 

Card  5  (315) 

ICATFL  =  Category  of  field 

=  1  for  material  property  random  field 
=  2  for  distributed  random  load  field 
=  3  for  random  structural  geometry* 

=  4  for  random  boundary  conditions* 

Note:  *  not  presently  available 

ITYPE  =  1  for  Young’s  modulus  random  field 

=  2  for  Poisson’s  ratio  random  field 
=  3  for  mass  density  random  field 
—  4  for  any  other  material  property  type 

ISORT  =  0  not  to  sort 

=  1  to  sort  eigenvalues 

(parameter  for  classifying  material  property  random  fields) 

Card  6  (15 ) 

ICDF  =  Marginal  cumulative  distribution  function  (CDF)  defining  the 

probabilistic  character  of  the  random  field 
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ICDF  Values  Corresponding  to  Various 
Probability  Distributions _ 


ICDF  Value 

Probability  Distribution 

1 

Fixed-valued  variable  j 

2 

Normal 

3 

One-sided  normal 

4 

Lognormal 

5 

Uniform 

6 

Triangular 

7 

Exponential 

8 

Beta 

9 

Chi-square 

10 

Gamma 

11 

Gumbel  (EVD-Type  I) 

12 

Maxwell  (not  presently 
available) 

13 

Rayleigh 

14 

Student-t 

15 

Weibull 

16 

Evd  Type  II  maximum 
(FrSchet) 

17 

EVD  Type  III  maximum 

Note:  EVD  =  Extreme  Value  Distribution 


Card  7  (2E10.3) 

RMV 

STD 

Card  8  (215,  E10.3) 
ICOREM 


Mean  value  of  the  field 
Standard  deviation  of  the  field 


Correlation  function  type  for  this  random  field 
0  for  fully  correlated  field 

1  for  exponential  correlation 

2  for  Gaussian  correlation 


P152713.PDF  [Page:  121  of  174] 


A 1-3 

=  3  for  triangular  correlation 

=  4  for  second-order  autoregessive  process  correlation 

IRHOPS  =  Parameter  that  defines  the  space  in  which  the  correlation  coefficient 

function  is  defined 

=  1  for  original  space  of  basic  random  variables  (x-space) 

—  2  for  transformed  standard  normal  space  (u-space) 

CORD  =  Correlation  distance  characterizing  the  random  field 

Card  8  (AS,  15) 

PREFX  =  Five  character  prefix  describing  the  SFEM  mesh  (i.e.  geometry  file)  to 

be  used  for  discretizing  this  field 

IPRF  =  0  no  physical  geometry  file  available  for  describing  this  SFEM  mesh 

=  1  PREFX. GOM  file  to  be  used  for  the  SFEM  mesh  (i.e.  identical  to 
mechanical  FE  mesh) 

=  2  PREFX.GOM  file  to  be  used  for  the  SFEM  mesh 

Card  9  (215) 

ICENT  =  0  for  automatic  computation  of  the  centroidal  coordinates  the  SFEMs 

=  1  for  user-supplied  centroidal  node  number 

=  2  for  user-supplied  global  coordinates  of  the  centroid  or  some  preferred 
location  for  discretization 

NSFEM  =  Number  of  stochastic  finite  elements 

_ _ Card  10  to  be  provided  NSFEM  times. 

Card  10a  (1615)  (provide  only  if  IPRF=2)  (  omit  if  ICENT=0  and  IPRF=1) 

JSFE  =  SFEM  number 

NFE  =  Number  of  mechanical  finite  elements  in  this  stochastic  finite  element 

JGRP  =  Element  group  number  of  the  mechanical  finite  elements  involved 

JFEM(k)  =  Element  numbers  for  the  (NFE)  mechanical  finite  elements 
(k=l,2, 

...,NFE) 

Card  10b  (1615)  (provide  only  if  IPRF=0) 


JSFE 


SFEM  number 
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A 1-4 

NNODES  =  Number  of  nodes  defining  the  boundary  of  the  stochastic  finite  element 

NBOUND(k)  =  Node  numbers  defining  the  NNODES  boundary  nodes 
(k=  1,2, 

NNODES) 

Card  11a  (15)  (provide  only  if  ICENT=1) 

INODE  =  Node  number  to  be  used  as  discretization  reference  point  (representing 

a  SFEM  centroid  or  some  other  preferred  location  within  the  SFEM) 

Card  lib  (3E10.3)  (provide  only  if  ICENT=2) 

XC(I)  =  Global  x-coordinate  of  discretization  point 

YC(I)  =  Global  y-coordinate  of  discretization  point 

ZC(I)  =  global  z-coordinate  of  discretization  point 

Note:  as  in  Card  1  la,  the  discretization  point  is  a  point  representing  a  SFEM  centroid 

or  some  other  preferred  location  within  the  SFEM. 
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APPENDIX  A 2 

INPUT  FORMAT  FOR  THE  SFEM-BASED  RELIABILITY 
ANALYSIS  PROGRAM:  STOVAST 
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INPUT  FORMAT  OF  PROGRAM  STOVAST:  PREFX.SFE 

Card  1  (A72) 

,  TITLE  =  Master  title  of  the  SFEM-based  reliability  analysis 

Card  2  (715) 

Seven  parameters  for  defining  the  master  control  codes: 

IRANFLD  =  0  No  random  field  definition  or  discretization  involved  in  this  analysis 

=  1  Random  field  discretization  involved  in  this  analysis 

IVBLDFN  =  0  No  definition  of  random  variables  required  other  than  those  generated 
from  random  fields 

=  1  Random  variables  are  defined  in  addition  to  those  generated  from 
random  fields 

IRESPNS  =  1  Limit  state  functions  (LSF)  or  response  functions  defined  in  terms  of 
VAST  finite  element  analysis 

=  2  LSF  defined  via  RESLIB  library 

ISYSDEF*  =  0  No  systems  defined  in  this  analysis 

=  1  Some  systems  are  defined  for  analysis 

ICOMREL  =  1  Component  reliability  analysis  via  FORM 

ISENSTY  =  0  Probabilistic  sensitivities  not  of  interest 

=  1  Computation  of  probabilistic  sensitivities  of  the  reliability  index  or 
failure  probability  are  of  interest 

ISYSREL*  =  0  No  systems  reliability  analysis  required 

=  1  Systems  reliability  analysis  required 

*Please  note  that  these  systems  are  not  yet  available  in  STOVAST.  Provisions  were  made  in  the  code 
to  permit  these  definitions  for  easy  extensions  of  STOVAST  capabilities  in  future  work. 


Omit  Card  3  if  IRANFLD  =  0 
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Card  3a  (A7) 


HEADR 

=  Seven  character  literal  constant  ‘IRANFLD’ 

Card  3b  (15) 

IOPT 

=  1  Random  variables  generated  from  random  fields  are  available  in 
PREFX.RVB  file 

_ _ _ _ _ _ _ _ _ 1 

Omit  Cards  4  and  5  if  IVBLDFN=0 

1 - - - 

Card  4a  (A7) 

HEADR 

=  Six  character  literal  constant  ‘IVBLDFN’ 

Card  4b  (15) 

NBRV 

=  No.  of  basic  random  variables  defined  for  this  analysis 

_  -  - - - - 

Provide  Card  5  NBRV  times) 

! - 

Card  5a  (A  8) 

VNAME 

=  Variable  name  (maximum  of  eight  characters) 

Card  5b  (215) 

ICDF 

=  Probability  distribution  function  for  this  random  variable 

(Admissible  values  of  ICDF  are  shown  below) 
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KRVT 


ICDF  Values  Corresponding  to  Various 
Probability  Distributions 

ICDF  Value 

Probability  Distribution 

1 

Fixed-valued  variable 

2 

Normal 

3 

One-sided  normal 

4 

Lognormal 

5 

Uniform 

6 

Triangular 

7 

Exponential 

8 

Beta 

9 

Chi-square 

10 

Gamma 

11 

Gumbel  (EVD-Type  I) 

12 

Maxwell  (not  presently 
available) 

13 

Rayleigh 

14 

Student-t 

15 

Weibull 

16 

Evd  Type  II  maximum 
(FrSchet) 

17 

EVD  Type  III  maximum 

Note:  EVD  =  Extreme  Value  Distribution 

=  0  Fixed  variable  (i.e.  deterministic  constant) 

=  i  FE  geometrical  nodal  coordinates  random  variable 
=  2  FE  material  property  random  variable  (eg.  E,  v,  p ) 
=  3  FE  size  property  random  variable  (eg.  A,  I,  t) 

=  4  FE  distributed  load  parameter  random  variable 
=  5  Fe  concentrated  load  parameter  random  variable 
=  6  Non-finite  element  threshold  random  variable 
=  7  Any  other  type  of  random  variable 
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Card  5c  (15) 

IOPX  =  1  Random  variables  defined  in  terms  of  mean  values  (/x),  standard 

deviations  (a),  and  bounds 

=  2  Random  variables  defined  in  terms  of  /x,  coefficient  of  variation  (COV), 
and  bounds 

=  3  Random  variables  defined  in  terms  of  distribution  parameters 

Card  5d  (3E10.3)  (required  only  if  KRVT^O) 

RMVL  =  Mean  value  (or  other  distribution  parameter) 

STDV  =  Standard  deviation 

BOUND  =  Bounds 

Card  5e  (E10.3)  (required  only  if  KRVT=0  or  ICDF=  1) 

FIXVAL  =  Fixed  value  of  the  variable 

Card  5f  (215) 

ICOREL  =  0  No  correlation  between  the  basic  random  variables 

=  1  Correlations  between  random  variables  defined  in  the  original  x-space 

=  2  Correlations  between  random  variables  defined  in  the  transformed 
standard  normal  u-space 

ICOPT  =  1  Correlation  described  in  terms  of  the  upper  triangular  portion  of  the 

correlation  matrix 

=  2  Correlation  defined  in  terms  of  correlation  pairs 
Card  5g  (omit  if  OCOREL=0  or  ICOPT=2)  (provide  J  times:  J  =  1,2,...,NBRV) 

Card  5g.l 

MM  =  No.  of  terms  in  row  number  J  of  the  correlation  matrix  starting  from  the 

diagonal  term  and  going  only  up  to  the  last  non-zero  term) 
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Card  5g.2  (free  format) 

AC(k)  =  Values  for  the  current  row  of  the  correlation  matrix 

(k=l,...,mm) 

Card  5h  (omit  if  OCOREL=0  or  ICOPT=  1) 

Card  5h.  1  (15) 

NCORPS  =  Number  of  correlation  pairs 

Card  5h.2  (free  format)  (provide  NCORPS  times) 

KI  =  First  variable  number  in  the  current  correlated  pair 

KJ  =  Second  variable  number  in  this  pair 

RHO  =  Correlation  coefficient  value 

Card  6a  (A7) 

HEADR  =  Seven  character  literal  constant  TRESPNS’ 

Card  6b  (15) 

IOPT  =  1  For  limit  state  (i.e.  component)  definition 

Card  6c  (15) 

NCOMPS  =  Total  number  of  components  to  be  defined  in  this  problem 

Card  6d  to  be  provided  NCOMPS  times. 

Card  6d  (A20) 

COMPNAME  =  Component  name  (maximum  of  20  characters) 

Card  6e  (415) 

NUMC  =  Identification  number  of  the  current  component 

LSFT  =  1  For  displacement-type  limit  state  function 

=  2  For  stress-type  limit  state  function 
=  3  For  buckling  limit  states  (not  presently  available) 
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NAVS 

Number  of  additional  random  variables  (i.e.  in  addition  to  those 
generated  from  random  fields)  to  be  used  for  defining  this  limit  state 

LTHR 

Number  of  the  random  variables  defining  the  limit  state  threshold 

Card  6f  (1615) 

II 

s-h; 

£  & 

Variable  numbers  of  the  NAVS  random  variables  for  the  limit  state 

NAVS) 

Card  6g  (215) 

LSDNN 

Limit  state  node  number 

LSDOF 

Limit  state  degree-of-freedom  number  for  displacements  or  stress  compo¬ 
nent  number  for  stresses 

Note:  The  convention  used  is  as  follows: 

Disnlacements 

Stresses 

u  =  1 

<7oc  =  1 

v  =  2 

0"yy  '  2 

w  =  3 

<Tzl  ~  3 

=  4 

<7xy  =  4 

<7y  =  5 

°yz  =  5 

ffz  =  6 

_  — -  "  '“1 

Omit  Card  7  for  ISYSREL=0. 

1 _  ■ 

Card  7a  (A7) 

HEADR 

Seven  character  literal  constant  TSYSDEF’ 

Card  7b  (15) 

NSYS 

Number  of  systems  defined  in  this  analysis 

_ _ _ _ - - -  1 

Provide  Card  7c  NSYS  times. 
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Card  7c.  1  (A20) 

SYSNAME  =  System  name  (maximum  of  20  characters) 

NCOMS  =  Number  of  components  defining  this  system 

Card  7c. 2  (215) 

ISYST  =  1  For  series  system 

=  2  For  parallel  system 

JSYS(k)  =  Component  numbers  of  the  components  that  define  this  system 

(k=l,  NCOMS) 

Card  8a  (A7) 

HEADR  =  Seven  character  literal  constant  ‘ICOMREL’ 

Card  8b  (15) 

IOPT  =  1  Component  reliability  analysis  based  on  the  Rackwitz-Ang-Tang  FORM 

algorithm  (not  presently  active  in  STOVAST) 

=  2  Component  reliability  analysis  based  on  the  HL-RF  FORM  algorithm 

=  3  Component  reliability  analysis  based  on  the  MHL-RF  FORM  algorithm 

Card  8b  (free  format) 

TOL  =  Tolerance  of  the  reliability  index  to  be  used  for  consequence 

ITERMAX  =  Maximum  number  of  iterations  permitted 

Card  8c  (free  format)  (omit  if  IOPT  =  l  or  2) 

CMERIT  =  Constant  c  in  the  definition  of  the  merit  function  (MHL-RF  algorithm). 

_ Omit  if  ISENSXTY =0. _ 

Card  9a  (A7) 

HEADR  =  Seven  character  literal  constant  ‘ISENSTY’ 
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Card  9b  (free  format) 


IPP 

=  1  Compute  parametric  sensitivities  of  jS  and  Pf  for  distribution  parameters 
only 

=  2  Compute  sensitivities  for  mean  values  and  standard  deviations 

=  3  Compute  1  and  2  above 

_ - _ — — - - - - - * - 1 

Omit  if  ISYSREL=0. 

L _  - 

Card  10a 

HEADR 

=  Seven  character  literal  constant  ‘ISYSREL’ 

Card  10b  (15) 

NSAN 

=  Number  of  systems  to  be  analyzed 

. .  _ _ _ _ _ — - “1 

Cards  10  to  lOe  are  to  be  provided  NSAN  times. 

- - — 

Card  10c  (1615) 

NSYS(I) 

=  System  numbers  of  selected  systems  to  be  analyzed 

(1=1,2, 

....NSAN) 

Card  lOd  (317) 

IB1 

=  1  Perform  systems  reliability  analysis  using  the  uni-modal  bound  technique 

IB2 

=  1  Compute  relaxed  bi-modal  bounds 

=  2  Compute  regular  bi-modal  bounds 

=  3  Compute  improved  bi-modal  bounds 

=  4  Compute  all  of  the  above 

IPNT 

=  i  Compute  PNET  estimation  of  the  failure  probability  for  this  system 

Card  lOe  (E10.3)  (omit  if  IPNET  ^  1) 

RHODEM  =  Demarcation  correlation  coefficient  to  be  used  in  the  PNET  procedure 
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(Note  that  systems  reliability  analysis  capabilities  are  currently  not  available  in  STOVAST.  So  the 
parameters  ISYSDEF  and  ISYSREL  should  be  set  equal  to  zero  for  all  STOVAST  analyses.) 
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STOVAST  OUTPUT  FILE  FOR  STOCHASTIC  BEAM  PROBLEM 


STOCHASTIC  FINITE  ELEMENT  BASED  RELIABILITY  ANALYSIS 
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TITLE  :  MAHADEVANS  STOCHASTIC  BEAM  PROBLEM 
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